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Chapter 1 


Introduction 


In 1985 Lagarias and Odlyzko [26] developed a general attack on knapsack cryp¬ 
tosystems which reduces solving subset sum problems to the problem of finding the 
Euclidean-norm shortest nonzero vector in a point lattice. Recent improvements to 
this attack [12, 19] have stimulated interest in finding lattice basis reduction algo¬ 
rithms well-suited to the lattices associated with subset sum problems. This thesis 
studies a new approach to lattice basis reduction originally developed by M. Seysen 
[38]. Seysen’s reduction algorithm was initially developed to find simultaneously good 
bases of a lattice and its dual lattice. However, it may also be successfully applied to 
solving subset sum problems, especially when combined with other known reduction 
methods. Using a collection of techniques, including Seysen’s algorithm, we show 
that it is possible to solve in practice a much larger class of subset sum problems 
than was previously possible. 


1.1 Point Lattices 

Let B be a set of vectors (bi, b 2 ,..., b n ) in R n . If these vectors are independent, 
then they form a basis of R n and any point x in n-space may be written as a linear 
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combination of vectors in B: 

n 

x = £ r,bi, for Tj e ]R, 1 < * < n. 

i—1 

Consider the set of points L C R n which may be written as the sum of integer 
multiples of the basis vectors: 

L = -I x = (aq,x 2 , ...,x n ): x = Yh z ibi, for z> € Z, 1 < * < n 
l «=i 

We call this set L the point lattice (or just lattice ) described by the basis B. 

Point lattices are pervasive structures in mathematics, and have been studied 
extensively. See [25], for example, for a survey of the field. In the area of combinato¬ 
rial mathematics alone it is possible to phrase many different problems as questions 
about lattices. Integer programming [20], factoring polynomials with rational co¬ 
efficients [27], integer relation finding [16], integer factoring [35], and diophantine 
approximation [36] are just a few of the areas where lattice problems arise. In some 
cases, such as integer programming existence problems, it is necessary to determine 
whether a convex body in R n contains a lattice point (for some specific lattice). In 
other cases the items of interest are short vectors in the lattice. As we shall see below, 
for certain cryptographic applications, we would like to be able to quickly determine 
the Euclidean-norm shortest nonzero vector in a lattice. 

It is important to note that the difficulty of finding the Euclidean-norm shortest 
nonzero vector in a lattice is an open question. If x = (aq,..., x n ), then the sup-norm 
of x, denoted ||x||oo, is defined as 


||x||oo = max |x;|. 

l<«<n 

It is known that finding the sup-norm shortest nonzero vector in a lattice is NP-hard 
[5]. Based on this evidence, we suspect that finding the Euclidean-norm shortest 
nonzero vector for any given lattice is also computationally difficult. However, it may 
be possible to find the shortest nonzero vector for many lattices quickly. Indeed, 
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current techniques for finding short vectors are slow in theory but often perform well 
in practice. 

The remainder of this chapter establishes the environment for our study of lattices 
and specific applications to cryptography. Section 1.2 discusses reduced bases of 
lattices and lattice reduction theory. Section 1.3 mentions some of the algorithms 
which currently exist for computing a reduced lattice basis B' given a basis B of 
a particular point lattice. In particular, we detail the operation of the Lenstra- 
Lenstra-Lovasz (LLL) basis reduction algorithm [27], which is currently the best 
known method for finding short vectors in a lattice. 

1.2 Reduced Lattice Bases 

Any lattice L may be described by many different lattice bases. Let Bi,B 2 ,... be 
distinct sets of vectors, all of which form bases of lattice L. We can imagine that 
there exists some ordering or ranking of the bases Bi, and thus one or more of the 
Bi might be considered “good” lattice bases of L. Lattice reduction theory deals with 
identifying “good” lattice bases for a particular lattice. If we are given a basis B 
which describes L, we would like to reduce B to basis B 1 , also describing L , where B 1 
is a “good” lattice basis in the sense of some reduction theory. 

There are two classical lattice reduction theories, one due to Korkin and Zolotarev 
[23, 24] and one due to Minkowski [29]. A basis B = (bi,b 2 ,... ,b„) of lattice L is 
said to be Minkowski-reduced if 

1. bj is the shortest nonzero vector in L. 

2. For 2 < ? < n, bi is the shortest vector in L such that (bj,...,bi) may be 
extended to a basis of L. 

Minkowski-reduced lattice bases always contain the shortest nonzero vector in the lat¬ 
tice. Subsequent basis vectors bi are selected by choosing the shortest lattice vector in 



4 


CHAPTER 1. INTRODUCTION 


L which is not a linear combination of the already selected basis vectors bi,..., bi_i. 
If bi = Y%= i ^jbj, Zj G Z, then it would be impossible to extend (bi,..., bi) to be a 
basis of L. 

The definition of Korkin-Zolotarev reduction is similar to that of Minkowski. We 
say a basis B = (bi,... ,b„) is Korkin-Zolotarev reduced if it satisfies the following 
three conditions: 

1. bi is the shortest nonzero vector in L. 

2. For 2 < * < n, let Si be the (i — l)-dimension subspace spanned by the basis 
(bi,... ,bi_i). Let Sj~ be the orthogonal complement of Si in R n . Finally, let 
Pi(L ) denote the orthogonal projection of L onto S^. Then choose bi such that 
P,-(bi) is the shortest nonzero vector in Pi(L). 

3. Size reduction condition. For 1 < i < j < n, 

|{C,(bi),/>.(bi))| < i||F i (bi)f. 

where Pi(x) = x. 

In the definition of Minkowski reduction, successive basis vectors bi are added to the 
lattice basis only if bi is the shortest vector in the lattice which will allow the basis 
to be extended. In Korkin-Zolotarev reduction, though, successive basis vectors bi 
are chosen based on their length in the orthogonal complement of the space spanned 
by the previous basis vectors bi,..., bi_i. 

Depending on the specific problem, we may find either, both, or neither of the 
above definitions of “good” lattice bases is sufficient. Certainly if the goal is to find 
the shortest nonzero vector in the lattice, as is the case with subset sum problems (see 
Chapter 3), we could use either Minkowski reduction or Korkin-Zolotarev reduction 
as a measure of how “good” a particular lattice basis is. If the item of interest involves 
multiple vectors in the lattice, one definition may be preferred over the other for that 
particular application. 
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1.3 Lattice Basis Reduction Algorithms 

Although both Minkowski reduction and Korkin-Zolotarev reduction provide frame¬ 
works for studying lattice basis reduction, computing a reduced lattice basis (in either 
sense) is in general a difficult problem. Currently, there are no known polynomial¬ 
time algorithms for finding either a Minkowski or a Korkin-Zolotarev reduced basis 
for a given lattice L. If such an algorithm existed, then we would be able to find 
the Euclidean-norm shortest nonzero vector in L in polynomial time by finding the 
reduced basis, for which b x is the desired vector. Thus, any polynomial-time lattice 
basis reduction algorithm we use will not be able to satisfy the strict conditions of 
Minkowski or Korkin-Zolotarev reduction. 

Techniques for finding relatively small vectors in a lattice have been known for 
some time (see [22] for example); it was not until recently, though, that a fast algo¬ 
rithm was known which was guaranteed to produce lattice bases with relatively short 
vectors. In [27] Lenstra, Lenstra and Lovasz described a polynomial-time algorithm 
for transforming a give lattice basis B = (b x \..., b„) of lattice L into an LLL-reduced 
lattice basis B' = (b x ,..., b[J. A basis B' is LLL-reduced if it has the following two 
properties: 


^ 2 

(b[, b?) 

for 1 < j < i < n, = ^ bf b ^, 

(1.1) 

br + > y b,*_ x 

for 1 < i < n, 

(1.2) 

where the parameter y G (^, 1) and b? = 

bj — 22l=i Mijb! (That is, B* = (b* x ,. 

• •5b*) 


is a Gram-Schmidt orthogonalized basis generated from B). Notice that LLL reduc¬ 
tion is similar to but weaker than Korkin-Zolotarev reduction. 


The Lenstra-Lenstra-Lovasz basis reduction algorithm converts lattice basis B 
into B' by performing two types of transformations. In a size-reduction step, LLL 
finds the largest j such that there exists an i > j and /qj violates Equation 1.1. By 
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performing the transformation 


bi < b; fJ bj? 

where |"-J denotes the nearest-integer function, we find that the new value of /qj is 
<5- In the exchange step, LLL searches for the smallest value of i such that bi_j 
and bi fail to satisfy Equation 1.2. Here LLL swaps the two vectors (bi_i <— b, 
and bi <— bi_i) to force compliance with the second LLL-reduced condition. The 
LLL basis reduction algorithm alternately performs size-reduction and exchange steps 
until both Equations 1.1 and 1.2 are satisfied, at which point the algorithm halts. 

LLL-reduced lattice bases satisfy a number of bounds. In particular, if y is the 
global LLL constant, then the first vector in the reduced lattice basis bi satisfies: 

||bi|| < (|^n) H X H’ for all x € L,x ± 0. 

In particular, for y = | (the value used in [27]), we have 

||t>i j| < 2 n_1 ||x||, for all x G L, x 0. 

Thus the length of bi is at most an exponential multiple of the length of the shortest 
nonzero vector in the lattice. (Similar bounds exist for the other vectors in the 
reduced basis.) In practice, the LLL algorithm usually performs much better than 
this exponential bound, although example lattice bases are known which cause the 
LLL algorithm to exhibit worst-case performance. 

Most of the work on lattice basis reduction algorithms since the introduction 
of LLL has focused on improving and extending the technique. For example, the 
version of LLL described above requires rational arithmetic (the variables in 
particular must be stored as fractions); multiprecision floating-point numbers are 
usually used to reduce the computation requirements, but may introduce error into 
the calculations. One method of reducing the multiprecision requirements is described 
in [33]. Similarly, [32] showed how to modify the LLL algorithm so that the set of 
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input vectors can be linearly dependent. Hierarchies of LLL-type algorithms have 
been investigated [34], stretching from LLL-reduction at one end of the spectrum to 
Korkin-Zolotarev reduction at the other. However, little effort has been expended on 
looking at algorithms not derived from or similar to LLL. 

This thesis examines an approach to lattice basis reduction of different structure 
than that of the LLL algorithm and related methods. This method was originally 
suggested by M. Seysen [38] to simultaneously produce a reduced basis and a reduced 
dual basis for some lattice L. Where the LLL algorithm concentrates on local opti¬ 
mizations to produce a reduced lattice, the Seysen approach considers the entire basis 
for methods of optimization. (Recall that the LLL size-reduction steps only consider 
Hi,j for the smallest possible j, and that only adjacent vectors bi_i and b; may be 
exchanged .) 

Chapter 2 describes the first phase of the research, in which Seysen’s basis reduc¬ 
tion algorithm and multiple variants were implemented and examined. Theoretical 
and empirical analyses of the fixed components of Seysen’s algorithm are given. For 
those parts of the algorithm which are permitted to vary, we examine some of the 
possible variations, and look at the effect of these changes on the performance of the 
algorithm. Possible extensions of the algorithm are also discussed. 

The motivation behind our study of Seysen’s lattice basis reduction algorithm is 
presented in Chapter 3. It is known that certain subset sum problems may be solved 
by finding the shortest nonzero vector in a particular lattice (the lattice is generated 
based on the specific construction of the subset sum problem). The best methods 
previously known for reducing subset sum problem lattices [26, 33] involve the LLL 
algorithm and some other heuristics, and are not very successful for n > 25 (n is the 
size of the subset sum problem to be solved). Chapter 3 details experiments which 
used Seysen’s algorithm in combination with the LLL algorithm and other heuristics 
to solve a much greater range of subset sum problems. 






Chapter 2 


The Seysen Basis Reduction 
Algorithm 


In 1990, Martin Seysen proposed a new method for performing lattice basis reduc¬ 
tion [38]. Seysen’s basis reduction algorithm (or just Seysen’s algorithm) differs from 
the LLL algorithm and its variants in that it considers all vectors in the lattice si¬ 
multaneously, and performs operations on those vectors which will reduce the lattice 
according to some measure. Recall that the LLL algorithm works locally on the lattice 
it is reducing; LLL will only perform an operation on two vectors which are adjacent 
to each other in the ordered lattice basis. 

Seysen was motivated to create a new method for basis reduction by a desire 
to find a better way to simultaneously reduce a lattice and its reciprocal (or dual) 
lattice. If lattice L is defined by basis vectors bi,..., b n , then the dual lattice L* of 
L is defined by basis vectors bj,..., b*, where 


( 2 . 1 ) 


(bi,bf) = l, 

(bi,bj) = 0, for i ^ j. 

Now consider what happens in the dual lattice when we perform a row move on bi 
and bj in L. (A row move is any operation which adds a constant multiple of one 
lattice basis vector to another basis vector.) Let bj = bj + Ab;, where A £ Z. We 
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consider what changes must occur in the dual lattice basis vectors bj,... , b*, since 
Equation 2.1 must hold at all times. For k ^ i, we find that: 

K' = bj, 

since 

(bj,b£) = (bj + Abi,b£), 

= (b j ,K) + A(b i ,b* k ), 

= 0 . 

For k = i, however, this is not the case: 

bf = bj* - Ab?, 

since 

(bJ,br) = (b j + Ab 1 ,bf / ), 

= (b j ,br / )+A(b i ,bf / ), 

= (bj,b? - AbJ) + A(bj, bj* - Abp, 

= (bj,K) - A(bj,bj) + A(bi,bD - A 2 (b,,bJ), 

= 0 — A + A — 0, 

= 0 . 

Thus, when we add Abi to bj in the basis of lattice L, we must subtract Abj from 
b- 1 in the basis of L*. It is easy to see that if lattice L is reduced with the LLL 
algorithm, the resulting reduced lattice may have a dual in which some basis vector 
is quite large, as no attempt is ever made to consider the size of the dual basis when 
row moves are being performed. Seysen’s algorithm attempts to choose row moves 
that reduce both the lattice and its dual. 
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We now outline the basic operation of Seysen’s algorithm. Let L and L* be a 
lattice and its dual which we wish to simultaneously reduce. Let A and A* be the 
associated quadratic forms of L and L*, respectively: 

^ = [ a 'j]l<t'J<n — [(k>’ kj)]l<ij<n , 

A * = K'Jl = [( bi,b j)]l<,',i<n' 

The element of matrix A is the inner product of the basis vectors bi and bj of 
lattice L. Notice that A and A* are inverses of each other and are symmetric. 

If L is the lattice defined by basis B, then any other basis B' of L may be written 
as: 

B' = BT, 

where T 6 SL n (Zi) (i.e. T is an n x n integer matrix with determinant 1). The 
quadratic form A' associated with B' may similarly be derived from A: 

A'= T* A T. 

For any quadratic form A, define the Seysen measure S(A) as follows: 

sw=Eo.x,=Eiib l riib:n j . (2.2) 

»=1 i=l 

A basis B is then S-reduced if: 

5(A) < 5(T* A T ), for all T € SL n ( Z). (2.3) 

We suspect that it is computationally difficult to find the optimal transformation 
matrix T for a given basis B. Consider, however, the class of transformation matrices 
T/\ defined by: 

- In A XUij, where i ± j, A G Z, 

I n is the n-dimensional identity matrix, and 

1 if k = i and l = j, 

Uij = K/]i< fc ,/<„ * where u ai = < 

0 if k yf i or / j. 
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(The matrix Uij has exactly one nonzero entry. Matrix TA has diagonal entries 
of 1 and exactly one nonzero off-diagonal entry.) Right-multiplying B by any TA 
simply adds A times the i th column of B to the j th column of B. If the columns 
of B are the basis vectors bi, then B TA is simply the transformed basis B' = 
(bx,b 2 ,... ,bj_i,bj + Abi,... ,b n ). 

Since it is easy to perform calculations with TA transformations, we focus our 
attention on products of one or more TA transformation matrices. It can be shown 
that every T £ SL n ( Z) may be written as such a product: 



We call a quadratic form S 2 -reduced if 


S(A) < S{TA A TA), for 1 < i,j < n, for A £ Z. 

Seysen suggests the following algorithm for ^-reducing a quadratic form: 
while (A is not S^-reduced) 
do 

choose i,j such that 3 A £ Z with 


5 (r,-i A T$ < S(A) 

let 


let 


A = 


\(°k 

2 Vli 


a = tAata 

where [-J denotes the nearest-integer function. This procedure for ^-reducing a 
quadratic form is Seysen’s basis reduction algorithm. 
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To date, little has been proven concerning the performance of Seysen’s algorithm. 
There are no known bounds on the Seysen measure of an .^-reduced basis (although 
bounds have been proven for .S'-reduced lattice bases), nor on the length of the short¬ 
est nonzero vector in the basis. The running time of Seysen’s algorithm is clearly 
bounded if the lattice basis consists of only integer vectors, but it is not known if the 
algorithm even terminates for basis vectors with real coefficients. However, prelimi¬ 
nary experiments performed by Seysen on lattices of dimension n < 30 suggest that 
this technique may be faster than the LLL algorithm and yield bases with shorter 
vectors. Based on these observations, a comprehensive investigation of theoretical 
and practical aspects of Seysen’s algorithm was undertaken. This chapter details the 
results of our study of Seysen’s basis reduction algorithm. 

Section 2.1 below discusses the theoretical underpinnings of Seysen’s algorithm. 
Empirical tests performed with various versions of the Seysen algorithm are detailed 
in Section 2.2. Section 2.3 mentions some modifications which may be made to 
Seysen’s algorithm when the performance of the basic version breaks down. Finally, 
Section 2.4 discusses possible extensions to Seysen’s algorithm. 


2.1 Theoretical Analysis 

We consider first the theoretical foundations of Seysen’s basis reduction algorithm. 
There are a number of questions concerning the actual structure of the algorithm 
which immediately arise. For a given quadratic form A, how might the ^-reduced 
and ^-reduced forms derived from A differ? Is it even sufficient to consider only TA 
transformation matrices, or are there lattices for which it is impossible to find the S'- 
reduced form using only TA transformations? How do we choose the order in which 
to apply the TA transformations, or equivalently how do we choose pairs of basis 
vectors for row moves? Is the Seysen measure function S'(A) = YJi=\ a i,i a *i,i a “good” 
way to rank different bases of a lattice? Finally, given that S{A) is an acceptable 
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measure function, is our choice of A = — “4)J ? given i and j, optimal? This 

section considers theoretical justifications for all of these questions. Section 2.2 below 
considers these questions from an empirical point of view. 


2.1.1 Sufficiency of Matrices 

As defined above, a basis B is S'-reduced if and only if its associated quadratic form 
A satisfies: 


S(A) < S(T f A r), for T G SL n { Z). (2.4) 


Thus, in order to Seysen-reduce a given lattice L which we know has basis B , we need 
to find a transformation matrix T G SL n ( Z) such that for all other T' G SL n ( Z) we 
have S(T* A T) < S{{T'Y A V). As SL n ( Z) is the set of all n x n matrices of unit 
determinant, we suspect that it is computationally difficult to find the desired matrix 
T directly. However, there are ways to avoid having to compute matrix T directly. 
Specifically, we can restrict our attention to a set of generating matrices for SL n ( Z), 
as we show below. 

Initially, let us assume that n = 2 and that A is the quadratic form associated 
with a lattice basis we wish to reduce. SL 2 CI 1 ) thus contains all 2 x 2 matrices [“ 
with ad — be = 1. Now, it is known [2, 37] that the set G 2 is a set of generating 
matrices for S'T 2 (Z), where: 


G2 — s 


1 A 
0 1 


agz 


U 


1 0 

A 1 


A G 


That is, if T G G 2 , then there exists a sequence of matrices 7\, T 2 ,..., T* such that 


T = Ti T 2 • • -T fc , T GG 2 for 1 < * < k. 
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(Actually, the set {[J }],[}?] : A G {-1,0,1}} is sufficient, since 


'l ±1_ 

X 

1 ±A 

[o 

l 


0 


h 

0 

X 

1 

«1 

±i l 


±A 1 


Section 2.2.2 below discusses the performance of Seysen’s algorithm when we restrict 


A = ±1.) 

Notice that the matrices [J *] and [} °] describe all possible row moves which can 
be performed on a 2 x 2 matrix. As an example, note that the matrix *So == [ — i o ] 
generated by: 


1 0 


1 1 


1 0 

-1 1 


0 1 


-1 1 


(So corresponds to swapping the two rows in a matrix.) Thus, the set of matrices 
TA for n = 2, i ^ j is exactly a set of generating matrices for SL 2 CZ). Therefore, it 
is sufficient for n = 2 for Seysen’s algorithm to consider only products of matrices of 
the form The difficulty is in choosing the right matrices and the right order of 
operations. 

Our analysis above assumed n = 2, but similar results are known where n is an 
arbitrary integer [30]. For fixed n > 0, 


G n = I n {J 


u 


1 


is a set of generating matrices for SX n (Z). Thus, it would be sufficient for Seysen’s 
algorithm to consider only TE and transformation matrices if it could pick the 
proper triples (i,j, ±1) at every step. In practice, Seysen’s algorithm chooses triples 
where A G Z, but the basic problem is still choosing the right triples in the 
right order. Choosing the correct (i,j) pairs to reduce and the value of A for that 
pair are discussed below. 
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2.1.2 Choosing Vector Pairs to Reduce 

Seysen’s algorithm does not specify how to choose which pair of basis vectors (b;, bj) 
to reduce on each iteration of the algorithm. At every iteration, it is necessary to find 
an (i,j) pair for which there exists a transformation matrix , A ^ 0, such that: 

S(T^AT t x tj )<S(A). 


Therefore, given that initially there are likely to be many pairs of vectors which may 
be reduced, we must decide how to select the best pair. 

Two options appear immediately as candidate vector selection methods: lazy 
selection and greedy selection. A lazy selection scheme simply chooses any available 
(i,j) pair in the easiest possible manner. For example, we can imagine two nested 
loops which generate (i,j) pairs and stop at the first pair for which A (i,j) =£ 0, where 


Kui) = 


M 

* 

3,3 


* 1,3 


*t,t / J 


Once such a pair is found, a transformation can be performed on the lattice 

basis. Then the algorithm could search for another ( i,j ) pair, perhaps continuing the 
search at the first pair lexicographically after ( i,j ). 

The second possible candidate selection method is a greedy approach. Here we 
calculate A(i, j, A) for each possible pair (i,j), where A(i,j,\) is defined: 


A(i,j,\) = S(T£,AT’: i )-S(A). 

Thus, any transformation matrix TA will have A(i,j, A) < 0. The algorithm then 
uses the pair of vectors (bi, bj) which minimizes A(i,j, A) in the next row move. 

One immediate disadvantage to a greedy approach is that it requires more exten¬ 
sive computations than the lazy selection method. This is true of any set of selection 
criteria which attempts to choose vector pairs to reduce in some fashion which per¬ 
forms better than random selection. If the two selection methods yield reduced bases 
of comparable Seysen measure, the added cost of an “intelligent” method may be 
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greater than the time saved by reducing the number of row operations. However, if 
one method should yield lattices with lower Seysen measure, the extra costs may be 
justified. 

We should point out that there is a distinction between choosing a pair of vectors to 
reduce and actually performing the reduction. Choosing a pair of vectors to reduce 
because they have the greatest potential to reduce the Seysen measure does not 
necessarily imply that we should perform the entire reduction and use the largest 
possible value of A. It may be wise to perform only a fraction of the possible row 
moves and reevaluate other possible pairs of vectors. We run the risk, if we are too 
greedy, of getting stuck too soon in a local minimum. 

There are reasons both to favor and to suspect the value of intelligent vector pair 
selection methods. One of the advantages that Seysen’s method has over the LLL fam¬ 
ily of basis reduction algorithms is that it looks at all the vector pairs simultaneously. 
The LLL algorithm works in a fashion similar to a bubble sort and LLL only considers 
row operations involving “adjacent” basis vectors (i.e. bi and bi_i for 2 < i < n). 
The cost of intelligent selection methods in terms of additional operations is certainly 
a disadvantage, but only if the cost is a significant fraction of the total running time. 
Section 2.2.1 below discusses these issues and presents empirical evidence of the cost 
and performance of a number of selection schemes for Seysen’s algorithm. From our 
experiments, the greedy selection scheme performs better than the lazy scheme, and 
the additional computation required to implement greedy selection is small. 

2.1.3 The S(A ) Function 

Equation 2.2 above introduced the Seysen measure function S(A) = 53”=i l|hi|| 2 ||bi l| 2 > 
the entire operation of Seysen’s lattice basis reduction algorithm centers on this quan¬ 
tization of relative reduction of two bases of lattice L. It is natural to question whether 
Seysen measures are indeed a reasonable means of ranking different lattice bases. We 
mention here some of the theoretical evidence which suggests that ranking lattice 
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bases by their Seysen measure is appropriate. 

The use of the quantity ||bi|| ||bjf || derives from elementary n-dimensional geometry. 
Recall the definition of the dual lattice B* = (bj,... , b*) of lattice B : 

(bi,bj) = 6 itj , for 1 < i,j < n, 

where Sij is the Dirac delta function (6,-j = 1 if i = j, Sij = 0 otherwise). Now, fix i 
and notice that 


(bijbj) = 1, 

l|bi || ||bj || cos(o;) = 1, 

llbillllbril = 


h, (2- 5 ) 

cos(a) 

where a is the angle between bi and bj. (Note that — ^ < oc < beacuse of the way 
in which b-‘ is defined.) 

Let Si denote the (n — l)-dimensional hyperplane spanned by the basis vectors 
bi,... ,bi_i,bi+i,... ,b n . Notice that b,* is perpendicular to Si by definition. Thus, 
given that a is the angle between b; and b-", the angle between bi and Si is 7r — ol. 
Thus, if ft — 7r — a, Equation 2.5 becomes 

INIIIKII = 

If basis vector bi is relatively dependent of the other vectors bj, 1 < j < n, j ^ *, 
then the angle between bi and the hyperplane Si will be relatively small, and thus 
will be large. Conversely, if b; is relatively independent of the other basis vectors, 


sin(/3) 


/3 will be close to | radians, and the product ||bi || Hb-* || will be close to one 1 . 

These geometric arguments lead directly to a measure which is a function of 
the products ||b; || ||bi || where 1 < i < n. Since \/3i\ < j, we could choose the 

as the measure function. Unfortunately, 


sin (Pi) 


function S'i(A) = E ||bi||llbf|| = E 

x Note that because of the duality between B and B *, we could also have considered /? to be the 
angle between b j and the hyperplane spanned by b£,..., bj*_ 1 , b; +1 ,... b* . 
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as Section 2.4.4 points out below, there is no simple formula for finding the optimal 
value of A for a row move involving the basis vectors bi and bj. Seysen is able to 
avoid these computational difficulties by using 

s(^) = ElNH|bni 2 , 

t=i 

n i 

_ S sin 2 (A)’ 

as the measure function, which does yield a simple formula for A. Since sin( ( 5 t ) £ [0,1], 
the squared terms in the S(A) function are guaranteed to be larger on a term-by- 
term basis than the corresponding terms in the Si (A) sum. Thus, if lattice basis B\ 
has smaller measure than basis B 2 using the Si measure function, B\ will also have 
smaller measure than B 2 when compared using the Seysen measure S. 

An additional advantage to using a function of the ||bi||||bj || product terms is that 
bounds exists on the size of the individual terms. In [17] Hastad and Lagarias show 
that the following bound applies for some primal basis B and dual basis B* of a given 
lattice: 

max { !|bi||, ||bj ||} < exp(0(n^)). (2.6) 

l<t<n 

This bound immediately implies that there exists a basis of L with Seysen measure 

1 

bounded by exp(0(n3)), since: 

max {||bi||, HKII} < exp(0(n*)), 

l<«<n 

max l|bi|ll|bi|| < exp(0(n5)) + exp(0(n3)) =exp(0(n^)), 

l<t<n 

l|bi||||bi || < nexp(0(n3)), 
i=i 

5(A) < exp(0(n3)). 

Seysen shows in [38] that the bound in Equation 2.6 may be improved to 


max {||bi||, ||bi*||} < exp(0((lnn) 2 )), 

l<«<n 


(2.7) 
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which reduces the bound on S'(A) for an ^-reduced lattice to: 

El|bi||||bf|| < raexp(0((lnn) 2 ), 

2 = 1 

< exp(lnn)exp(0((lnn) 2 ), 

S(A) < exp(0((lnn) 2 ). 

To date, this is the best known bound on the Seysen-measure of an ^-reduced lattice 
basis. However, as is the case with the LLL algorithm, in some cases Seysen’s algo¬ 
rithm produces ^-reduced lattice bases which have measures much lower than the 
theoretical bound. 


2.1.4 Choosing A Values 


We now consider the choice of A values in Seysen’s basis reduction algorithm. Assume 
that S'(A) is as in Equation 2.3 above, and that only two-vector row moves are 
considered (i.e. transformation matrices of the form T— for integer values of i,j and 
A). We first show that 


1 




2 




3,3 



Cli,i 


( 2 . 8 ) 


yields the maximum possible reduction in S(A) for fixed values of i and j where 
A 6 R. Further, we show that if we require A E Z, then 


1 ( ai 'i — Cli, i 

“ \ a *j,j a M 


(2.9) 


is indeed the value which yields the best possible reduction in the Seysen measure of 
the lattice basis. 

Let i,j be fixed integers with 1 < i,j < n; bi,bj,bj* and bj are the basis vectors 
on which we will perform a row move. Without loss of generality, we assume that Abi 
will be added to bj and Abj will be subtracted from b-\ Define bj and b^ 7 to be the 
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values of bj and b,* after the row move is performed: 

bj = bj + Ab;, 
bi / = b i *-Ab j *. 

Let A and A* be the quadratic forms associated with the lattice and its dual before 
the row move occurs, and let A' and A*' be the associated quadratic forms after the 
row move. Then 


A! = Tti A T x 


A*' = T~ x A* T~ a . 


Now, given that T x - transition matrices have exactly one off-diagonal nonzero 
entry, it is easy to see that A' differs from A only in the values in the d row, the 
j th row, the f th column, and the j th column. The same is also true for A* 1 . Since the 
Seysen measure function S(A) only depends upon the diagonal elements in A and A , 
we know that 




i =1 


S = 1 


- a i,i a i,i A a j,j a j,j a 3,3 a j,j * 


( 2 . 10 ) 

( 2 . 11 ) 


When A is right-multiplied by T) A -, the i th column of A is added to the j th column of 
A. When this matrix is subsequently left-multiplied by TA, the row is added to 
the j th row. Thus, after these two transformations, the value of a,-,,- is unchanged, but 


3,3 


= a 


3,3 


+ 2Aa,y T A 


( 2 . 12 ) 


If we perform a similar analysis in the dual quadratic form A*, we find that 
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Using Equations 2.12 and 2.13, Equation 2.11 becomes: 

5(A') - 5(A) = a iti (<,- - 2A < j A A 2 a*+ aY (a jtj + 2Xa i>j + A 2 a ijt ) 

* * 

= 2A A 2A aijdjj 2\di t id it j. 


Differentiating with respect to A and setting the result equal to 0, we find that: 


^ (S(A') - 5(A)) = 4A a iti aY A 2 a^a*- - 2a t> < J = 0, 

4Aoj i (fl,- j = 2di t {d-j 2ot i,jdjjf 


A = 


a i,i a i,j a i,j a j,j 


2 diid* ,■ 


3,3 


1 

2 


<d 


3,3 


**,3 




Thus, if A could take on any real value, for fixed i and j the minimum value of 5(A') 


is obtained with A = | 


aj,; 
a, , 


We have shown that the minimum value of 5(A) with A G R is obtained when 
A satisfies Equation 2.8. Our goal now is to show that if A is restricted to integer 
values, Equation 2.9 yields that value of A for which 5(A) is minimized. Let 


A = X z A A r , where A* € Z,0 < A r < 1, 

A(i,j,\) = S(T* i AT? J )-S(A). 

We know that for fixed t,j,A = ± (A - f“) minimizes the value of A. Furthermore, 
as A is a quadratic function of A, at least one of X z and \ z + 1 must minimize A for 
fixed, integer values of A. 

Consider A(i,j,X z ) and A(i,j,X z A 1): 

A(^,j, A^) — 2A z di^djj T 2\ z dijd2\ z di^d^j) 

A(i,j, A* + 1) = 2(X Z + l) 2 di t id* j A 2(A Z + 1 ) a i,j a j,j — 2{X Z A l)^*,***,-^, 

= 2A z d{ t idjj A 2X z d{jdjj 2X z di'id { j 

A lX z di t {dj j A 2(ii.i A 2ctt, j&jj Set; 
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Thus, 


A(i, j, \ z + 1) — A(*, j, A*) — 4A z ai t ia*- - + 2Oi^jj + 2a,ja J j 2a,- i ;a i> j. (2.14) 

As S( A) is a non-negative valued function which we want to minimize, we are interest¬ 
ed in large, negative A values (i.e. |A| should be large, A < 0). Thus, if Equation 2.14 
is > 0, we should choose A = X z ; similarly, if Equation 2.14 is < 0, set A = A* + 1. 
When is Equation 2.14 greater than zero? 


A(i, j, A z + 1) — A (i,j, X z ) > 0, 

4:X Z T 2u^icijj T 2ciij(ij j 2o,i t iCi- j !> 0, 

4A z o,i t iQ,j j > 2d{i(i^j 2o.i t i<ij j 2(ii jQj j , 

- 1 ), 

2 \ a h a ij J 
x z > x -§, 

X z X z -(- x r 2 > 

-^ A r < 2- 


Thus, if A r < 1, then Equation 2.14 is positive, and we should choose A' = X z . If 
A > Equation 2.14 is negative, and we should set A' = X z + 1. Thus, 


X' — \X Z + A r J , 

— i ( Q *’ J _ ai 'i 

2 \ a *hi Qi ’ i 


which proves Equation 2.9. 


2.2 Empirical Analysis 

In the previous section we attempted to provide some theoretical justification for 
Seysen’s basis reduction algorithm. The basic analysis suggests that Seysen’s tech¬ 
nique is viable, but as yet there are no significant bounds on the running time of 
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the algorithm. In this section we detail numerical experiments which were performed 
using Seysen’s algorithm. These experiments yield greater insight into the classes of 
lattices best suited for reduction by Seysen’s algorithm, as well as an indication of 
the effectiveness of Seysen’s technique. 

Before we begin detailing the empirical results it is appropriate to detail our 
test conditions. All implementations of Seysen’s basis reduction algorithm and the 
LLL algorithm were written in FORTRAN. Multiprecision floating point arithmetic 
operations were performed by a package of routines written by David Bailey at the 
NASA Ames Research Center [4]. Tests were run on Silicon Graphics, Inc., IRIS-4D 
workstations; the IRIS uses the MIPS R3000 chip set as its main processor. 

The experiments described below explore many of the same aspects of Seysen s 
algorithm discussed in the previous section. Section 2.2.1 compares lazy and greedy 
schemes for choosing the row move to perform on each iteration of the algorithm. The 
effects of restricting A choices are discussed briefly in Section 2.2.2. Sections 2.2.3 
and 2.2.4 compare the performance of Seysen’s algorithm with the LLL lattice on two 
classes of lattice bases. 


2.2.1 Lazy vs. Greedy Selection Methods 

In Section 2.1.2 above we outlined the differences between lazy and greedy methods 
for selecting a row move to perform on each iteration of the Seysen while loop. In this 
section we describe the results of test lattice reductions on small-dimension lattices 
using both the lazy and greedy approaches. 

All experimental tests were performed on random integer lattices of determinant 
one. Integer lattice bases were generated as follows. Let B be an n x n integer matrix 
where the i th column of B corresponds to basis vector bi. The elements of matrix B 
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Table 2.1: Comparison of Lazy and Greedy Selection Methods 


n 

Avg. # Steps (Lazy) 

Avg. # Steps (Greedy) 

Ratio (Lazy/Greedy) 

20 

2079.90 

758.50 

2.74 

25 

4096.40 

1624.25 

2.52 

30 

7444.80 

3279.45 

2.27 

35 

8787.35 

3094.25 

2.84 


are: 


1 


if i = j, 


B = [M 




= < 


rand(x) — \x 


if * >3, 

if i < j- 


The function rand(x) generates a random number chosen uniformly from the interval 
[0, x]; in these experiments x = 4. Notice that det(i?) = 1 since B is upper-triangular 
and all diagonal entries & ti ; = 1. 

To generate a random, dense lattice which is not upper-triangular yet still has 
determinant equal to 1, we perform random row moves on matrix B to generate 
matrix B'. We choose n 2 pairs of integers (i,j) with 1 < i,j < n and i ^ j. For each 
such pair, A is chosen to be +1 or —1 with equal probability. Then, the current b; 
is scaled by A and added to bj. (That is, we set B = B T&.) The result of these 
n 2 random row moves is matrix B 1 , which is a basis for our test lattice L. Since TA 
transformations preserve the determinant of the lattice, we know that det(L) = 1. 
We may thus measure the performance of Seysen’s algorithm on lattice L by how 
close reduced lattice L' is to I n . 

Table 2.1 summarizes the results of tests comparing the performance of lazy and 
greedy selection methods. Twenty test lattices were generated for each dimension 
n G {20,25,30,35}. In all cases where n < 30, both lazy and greedy algorithms were 
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able to completely reduce all test lattices to I n . The table shows the average number 
of row moves required by the lazy and greedy methods to reduce a lattice to I n . On 
average, the lazy selection scheme required over twice as many row reductions as the 
greedy scheme did to reduce a given test lattice. 

At n = 35 the lazy algorithm was able to reduce to I n only two of the twenty 
attempted lattices; the remaining problems all encountered local minima during the 
reduction, thus halting Seysen’s algorithm. The greedy implementation was unable to 
completely reduce any of the n = 35 test lattices. The two versions of the algorithm 
performed about equally well if we look at the Seysen measure of the reduced n = 35 
test lattices. 

These experimental results tell us two things concerning the relative merits of lazy 
and greedy selection schemes. First, when both lazy and greedy methods are likely 
to produce lattice bases with similar Seysen measures, the greedy selection methods 
will save at least a factor of two in the number of reduction steps. Second, based on 
the n = 35 data, using greedy instead of lazy does not appear to significantly reduce 
the performance of the algorithm as a whole. For our test lattices neither method 
performed significantly better than the other in terms of the Seysen measure of the 
^-reduced lattice bases. 

One might argue that it is not reasonable to compare only the number of reduction 
steps required to reduce lattices using greedy and lazy selection methods, since that 
measure fails to take into account the cost of selecting the two basis vectors to reduce. 
A naive implementation of the greedy algorithm might require 0(n 2 ) time, as there 
are \n(n — 1) possible pairs of basis vectors (bpbj),*' ^ j which must be considered. 
However, it turns out that, after an initial 0(n 2 ) precomputation phase, only O(n) 
time is required to greedily select the next row move. Assume that we have computed 
A(i,j,X(i,j)) values for all pairs of integer (z, j), 1 < i,j < n. Now, after a specific 
row move involving basis vectors i' and j' is performed, the only previously computed 
values of A which need to be updated are those for which i = i',i = j',j = i' or 
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j = j'. (If you consider A to be an array of values, the and (j ) rows and 

columns of A are all that need to be recomputed.) Thus, this recomputation can be 
performed in 0(n ) time. 

Storing A values can reduce the cost of a greedy selection method from 0(n ) to 
0(n). However, even 0{n) cost would be prohibitive if the actual amount of compu¬ 
tation required to select a pair of vectors was comparable to the cost of performing 
a row move. This is not the case; performing a row move requires 0(n) multipreci¬ 
sion math operations, whereas the stored A values need only be stored as single- or 
double-precision floating point numbers. (The values are generally different enough 
that 32 or 64 bits will provide more than enough precision.) Since multiprecision op¬ 
erations take significantly more time than double-precision operations, and since each 
row-move requires 0{n) operations, it seems reasonable to discount the added cost of 
performing the greedy selection as noise when compared to the cost of implementing 
the main portion of the algorithm. 

Based on these experimental results, we chose to use a greedy selection strategy 
in all subsequent implementations of Seysen’s basis reduction algorithm. For lattices 
where we know that Seysen’s algorithm will be able to perform a significant amount 
of basis reduction, such as the sparse lattices associated with subset sum problems 
(see Chapter 3), the greedy selection method and its expected reduced execution time 
are preferred. 

2.2.2 Choosing A Values 

As shown in the previous section, the selection method for choosing row moves in 
Seysen’s algorithm can affect both the algorithm’s performance and running time. In 
our comparison of lazy and greedy selection methods above, however, we implicitly 
sidestepped another such issue: the method by which values of A are chosen for a 
specific row move. Section 2.1.4 showed that for specified lattice basis vectors (bp bj), 
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the value of A which minimizes A (i,j, A) is: 



T 

K, 

a i,j \ 


2 


a i,i ). 


(2.15) 


However, recall from Section 2.1.1 that any transformation matrix T-j may be repre¬ 
sented as the product of A Tg 1 matrices (+1 if A > 0, -1 otherwise). Thus, when a 
T a - transformation is applied to lattice basis B, we are implicitly performing A row 
moves at once. It may be the case that for some lattices, a finer degree of control is 
required; that is, a greedy algorithm might perform even better if it was restricted to 
performing only Tf j and T~J transformations. That way the algorithm would have 
the finest possible degree of control over row operations. 

It is also important to note that a greedy selection mechanism uses A values in 
two distinct places. First, in order to select a pair of basis vectors to reduce, the 
greedy approach calculates A.(i, j, \{i, j)) for all possible values of i,j,i ^ j (A(z, j ) 
is the function in Equation 2.15 above). Once a pair of vectors has been chosen, 
a T— transformation is applied. In the first case, A is used as part of the scoring 
mechanism in order to choose a set of basis vectors to reduce. In the second case A 
plays a different role, the number of times to add vector bj to bj. Because A values 
have these two distinct functions, it is important that we distinguish between those 
roles when testing methods of choosing A values. 

We consider in this section three versions of Seysen’s basis reduction algorithm 
and their performance on a set of randomly generated integer lattices 2 . All three 
versions use a greedy selection scheme to choose the next row move to perform. They 
differ only in the set of allowed values of A in the scoring and transformation phases. 
These version are: 


1. (Z, Z): A may take on any integer value when choosing a set of basis vectors for 
the next row move, and any T)j may be performed on those vectors. 

2 In fact, we use the same set of integer lattices used in the previous section for comparing the 
lazy and greedy selection mechanism. 
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2. (Z, ±1): A may take on any integer value when choosing a set of basis vectors 
for the next row move. However, only T}- and T~J actual transformations are 
allows. (If A > 0 we add b; to bj. If A < 0 we subtract b; from bj.) 

3. (±1,±1): A may only be ±1 when choosing the next row move, and only T^ 1 
transformations may be performed on the basis. 

The (Z, Z) version is identical to our “greedy” implementation of the previous section; 
it will serve as our control. The (Z,±l) version of Seysen’s algorithm is designed to 
greedily select the best possible row move based on unlimited A values, but to perform 
the least possible number of changes to the lattice basis before recomputing what to 
do next. The (±1,±1) version also restricts lattice basis changes to the minimum 
amount possible on each step, but this version selects a row move based only on what 
it can do immediately to reduce the S(A) measure, not on any “future potential.” 

Table 2.2 compares the performance of the (Z, Z), (Z, ±1) and (±1, ±1) versions of 
Seysen’s basis reduction algorithm. For each value of n, twenty test lattices of dimen¬ 
sion n were generated and Seysen-reduced by each of the three methods. The table 
lists aggregate information for each value of n (all numbers shown are the geometric 
mean of the experimental values obtained for each of the twenty test lattices): 

• n : The dimension of the test lattice in this group. 

• L\q\ E"=i l°gio IIbill before reduction. 

• -S'(H): The Seysen measures before reduction. 

• L* 10 : E?=i logician before reduction. 

• Method: The restructions placed of A values. 

• L' w : E?=i l°giol|bi|| after reduction. 


• The Seysen measure after reduction. 
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Table 2.2: Comparison of (Z,Z), (Z, ±1) and (±1,±1) Options 


n 

Lio 

S(A) 

r* 

■^10 

Method 

T> 

-^10 

S(A') 

r * / 
^10 

# Steps 





(Z,Z) 

0. 

20. 

0. 

757.4 


82.4 

4.7 • 10 11 

118.2 

(Z,±1) 

0. 

20. 

0. 

767.0 





(±1,±1) 

0. 

20. 

0. 

787.3 





(Z,Z) 

0. 

25. 

0. 

1622.1 


132.0 

5.7 • 10 14 

192.9 

(Z,±1) 

0. 

25. 

0. 

1603.4 





(±1,±1) 

0. 

25. 

0. 

1610.3 

m 




(Z,Z) 

0. 

30. 

0. 

3275.3 

30 

195.8 

9.3 • 10 17 

287.2 

(Z,±l) 

0. 

30. 

0. 

3176.1 





(±1,±1) 

0. 

30. 

0. 

3316.1 





(Z,Z) 

110.4 

1.0 • 10 8 

112.6 

3037.0 

35 

269.2 

6.2 • 10 20 

390.0 

(Z,±l) 

99.0 

4.0 • 10 7 

102.2 

3022.0 





(±1,±1) 

112.3 

1.3 - 10 s 

114.6 

2920.8 





(Z,Z) 

216.9 

5.2 • 10 12 

227.1 

2240.2 

40 

343.8 

1.8 • 10 23 

507.0 

(Z,±l) 

206.5 

1.6- 10 12 

217.0 

2399.4 





(±1,±1) 

205.0 

1.5 -10 12 

217.3 

2527.2 





(Z,Z) 

304.6 

4.2- 10 15 

323.5 

2369.9 

45 

434.0 

2.0 • 10 26 

656.7 

(Z,±l) 

290.6 

1.2 • 10 15 

312.3 

2569.8 





(±1,±1) 

296.5 

1.9 • 10 15 

315.8 

2739.2 





(Z,Z) 

402.8 

2.3 • 10 18 

429.2 

2698.1 

50 

541.2 

4.5 • 10 29 

833.5 

(Z,±l) 

386.1 

6.3 • 10 17 

418.4 

3276.4 





(±1,±1) 

393.8 

1.1 • 10 18 

422.4 

3491.1 
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• L\ o= E"=il°gioll b ill after reduction. 

• # Steps: The number of row moves performed during the reduction. 

For n < 30, all three implementation were able to completely reduce all test lattices 
to I n . The only difference in the performance of the three methods was in the number 
of reduction steps required to reduce a test lattice, and these differences were minor 
(no more than 5% variation among any of the values for a given dimension). 

More differences among the methods appeared once n > 35 and Seysen s algorithm 
was no longer able to reduce any of the test lattices to I n . For the majority of the test 
lattices, the (Z, ±1) appears to yield the most Seysen-reduced lattice basis, although 
it requires significantly more row moves to perform this reduction than the (Z,Z) 
method. This improvement is expected; after all, the only difference between the 
(Z, Z) and (Z, ±1) methods is that that latter looks more frequently at the lattice basis 
and takes smaller “steps” while performing a reduction. However, as the dimension 
of the lattice basis increased, the ratio of row moves required by (Z,±l) to row 
moves required by (Z,Z) also increased. By the time n = 50, the (Z,il) method 
required approximately 30% more reduction steps to reduce test lattices than the 
(Z,Z) method did. 

The performance of the (±1, il) method fell somewhere between that of (Z, ±1) 
and (Z, Z) for test lattices with n > 45. This is somewhat surprising in and of 
itself, since the capability of the (±1,±1) method to consider future reduction is 
severely limited. This unexpected performance may be an artifact of our method of 
generating test lattices; we only performed n 2 row operations on the initial upper- 
triangular lattice basis and used only transitions to modify the lattice basis. This 
could account for the relatively good performance of the (±1,±1) method, and also 
the differences between the (Z,Z) and (Z,±l) methods. 

The experiments summarized in Table 2.2 do not indicate that any one of the 
tested methods consistently performs better than the others. Without any clear 
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indication that one method would yield significantly better results (especially as n -> 
100), we were reluctant to use any method in Seysen’s algorithm except the (Z, Z) 
one implicitly defined by Seysen himself. For special types of lattices, it is probably 
worthwhile to compare these methods again. It may be that increased performance 
by the (Z, ±1) method more than offsets the increase in the number of row operations. 
However, for our experiments in subsequent sections (and the subset sum lattices of 
Chapter 3) we continued to use the (Z,Z) form of Seysen’s algorithm. 

2.2.3 Testing the Be Lattice 

The previous sections detailed small tests which compared the relative performance of 
Seysen’s algorithm when a few minor changes were made to its structure. In this sec¬ 
tion and the following one we investigate more fully the performance of the algorithm 
itself as a reduction technique for self-dual lattice bases and random integer lattices. 
Our next series of tests was suggested by J. C. Lagarias (personal communication) to 
test the self-dual reduction performance of the algorithm. Let B$,0 < 0 < | be the 
lattice basis consisting of the following basis vectors: 

b 1 = (l,0,0,0,...), 
b 2 = (0,l,O,O,...), 
b 3 = (-M,i,0,...), 

b 4 = (0,-M,!,•••), 


If B$ is represented as a matrix with b; as the i th column, then we have: 


B, 


1 


if i = j, 


[hi\ 




< 0 if i > jy 

(—1 y~ t+1 0 if i<j. 
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Basis B$ has the property that ||b-‘|| grows exponentially with n, but rather slowly 
for small dimensions. 

Tests were performed on Bq lattices with 0 = 0.4 using Seysen’s basis reduction 
algorithm for dimensions 5 < n < 105. Based on the experimental results of Sec¬ 
tions 2.2.1 and 2.2.2 above, we used a greedy selection method to choose which pairs 
of vectors to reduce on each iteration of the algorithm, and A was allowed to be any 
integer value during all stages of the algorithm. The results of these tests are given 
in Table 2.3. For each test lattice, the following information is given: 

• n: The dimension of the lattice. 

• Lio : E?=il°giol|bi|| before reduction. 

• The Seysen measure of Bg before reduction. 

• L* 10 : £”=i logroll KII before reduction. 

• L' 10 : E"=i log 10 11bi|| after reduction. 

• 5(^4'): The Seysen measure of Bg after reduction. 

• Tio : £"=i l°gio ||b||| after reduction. 

• # Steps: The number of row moves performed during the reduction. 

The exponential growth of ||bj , || may be easily seen by looking at the rate of growth of 
the L\ q column in the tables. Remember that L\ 0 is the sum of the base 10 logarithms 
of the lengths of the dual basis vectors. This sum grows at least linearly with respect 
to n; thus, the ||bj || grow exponentially in n. 

For the Bg lattice Seysen’s basis reduction algorithm yields little improvement 
in the vector lengths ||bi||. Indeed, for some values of n we have Lio < L' w . This 
is not the case, though, for the dual lattice; Seysen’s reduction algorithm greatly 
decreases the lengths of the vectors in the dual basis. When the algorithm completes, 



34 


CHAPTER 2. THE SEYSEN BASIS REDUCTION ALGORITHM 


Table 2J 


Performance of Seysen’s Algorithm on B$ for 6 — 0.4,5 < n < 105 


n 

Lio 

5(A) 

r* 

-^10 

T> 

-^10 

S(A') 

r* / 
•^10 

# Steps 

5 

0.30 

7.36e+00 

0.50 

0.30 

6.59e+00 

0.30 

3 

10 

1.10 

3.65e+01 

3.80 

1.10 

1.61e+01 

1.00 

12 

15 

2.30 

1.88e+02 

10.60 

2.10 

2.78e+01 

1.90 

29 

20 

3.70 

1.00e+03 

21.10 

3.40 

4.28e+01 

3.10 

45 

25 

5.30 

5.37e+03 

35.20 

4.90 

5.80e+01 

4.20 

70 

30 

7.10 

2.89e+04 

53.00 

6.30 

7.88e+01 

6.20 

135 

35 

9.10 

1.55e+05 

74.40 

8.20 

1.03e+02 

8.00 

177 

40 

11.20 

8.35e+05 

99.50 

10.40 

1.30e+02 

9.90 

240 

45 

13.40 

4.49e+06 

128.30 

13.10 

1.72e+02 

12.80 

363 

50 

15.70 

2.42e+07 

160.70 

15.80 

2.32e+02 

16.90 

509 

55 

18.20 

1.30e+08 

196.70 

18.60 

2.62e+02 

18.30 

698 

60 

20.70 

6.99e+08 

236.40 

22.30 

3.53e+02 

22.50 

772 

65 

23.30 

3.76e+09 

279.70 

25.00 

4.01e+02 

25.70 

1067 

70 

25.90 

2.02e+10 

326.80 

26.40 

4.40e+02 

28.10 

1467 

75 

28.70 

1.09e+ll 

377.40 

32.70 

6.93e+02 

36.60 

1702 

80 

31.50 

5.85e+ll 

431.70 

31.90 

6.00e+02 

35.70 

2123 

85 

34.40 

3.15e+12 

489.70 

35.70 

7.25e+02 

40.40 

2775 

90 

37.30 

1.69e+13 

551.30 

39.70 

8.50e+02 

45.40 

3345 

95 

40.30 

9.10e+13 

616.60 

44.10 

1.03e+03 

51.00 

4839 

100 

43.30 

4.89e+14 

685.60 

49.00 

1.20e+03 

55.50 

6363 

105 

46.40 

2.63e+15 

758.10 

50.80 

1.19e+03 

56.70 

8303 



2.2. EMPIRICAL ANALYSIS 


35 


700t 


SOOf 


500 


10 





400 




36 


CHAPTER 2. THE SEYSEN BASIS REDUCTION ALGORITHM 


we have L' 10 « L* 0 ' and the lengths of the vectors in the prime and dual lattice bases 
are comparable. Figure 2-1 shows graphically the improvement made by Seysen’s 
algorithm to ||bj ||. 

The results obtained from application of Seysen’s algorithm to Bg lattices are 
quite promising. The algorithm was able to significantly reduce the lengths of the 
dual basis vectors b- 1 without significantly increasing the lengths of the basis vectors 
for Bg themselves. In fact, the resulting primal and dual bases have basis vector 
lengths which are comparable. Certainly this suggests that Seysen’s algorithm is a 
viable technique for applications in which we wish to simultaneously reduce a lattice 
and its dual. 


2.2.4 Testing Random Integer Lattices 

The reductions of Bg lattices tested application of Seysen’s algorithm to a very narrow 
class of lattices. From a cryptographic point of view (see Chapter 3 below), and in 
many other cases, our goal is to reduce random lattices with integer basis vectors. 
In general, we do not know that our lattice conforms to some specific structure; we 
are give only the basis vectors themselves. Thus, it is appropriate to investigate the 
performance of Seysen’s algorithm on randomly generated integer lattices. 

When we considered in Section 2.2.1 the choice of a lazy or greedy selection method 
for Seysen’s algorithm, we ran tests on random integral lattices L with det(Z) = 1 
of small dimension (n < 35). We utilize the same technique for generating random 
lattices as used before, except now we consider lattices up to dimension 50. In addition 
to running Seysen’s algorithm over these test cases, each lattice was also reduced using 
the LLL algorithm with y < 1.0. 

Table 2.4 summarizes the results of these experiments. For each value of n = 
0 mod 5, 20 < n < 50, twenty different test lattices were generated. The columns n, 
iio, S(A ), LI o, L[ o, S'(A'), and L\ 0 ' identify the same quantities as in Table 2.2 above. 
The column labeled Steps (Seysen)” reports the average number of reduction 
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steps (row moves) performed by Seysen’s algorithm for each test lattice. The average 
number of row operations performed by the LLL algorithm is listed in the Steps 
(Lovasz)” column. (LLL reduction was not performed on lattices with n > 45 because 
of the excessive amount of computer time required to obtain results). 

The LLL basis reduction algorithm was able to reduce all tested lattice bases to 
the n-dimensional cubic lattice basis (i.e. B' = J n ), which has Seysen measure zero. 
Seysen’s algorithm performed similarly on all lattice bases with n < 31. (Values in 
parentheses in the Z/ 10 column indicate the number of lattice bases which were Seysen- 
reduced to the n-dimensional cubic lattice.) For 32 < n < 34 the Seysen algorithm 
was able to completely reduce only some of the attempted lattice bases; no lattice 
basis with n < 35 was ever reduced by the Seysen algorithm to I n . 

The degradation in the performance of Seysen’s algorithm as n increases from 
30 to 35 is quite surprising. Apparently, for these types of lattices, the probability 
of reaching a lattice basis which is a local minimum of the Seysen measure function 
increases substantially over that range. The LLL reduction algorithm does not exhibit 
any decrease in performance, except for an overall increase in the number of reduction 
steps required to convert the given lattice basis into the n-dimensional cubic lattice 
basis. Of course, the LLL algorithm does take significantly longer to run than the 
Seysen algorithm, but these tests suggest that Seysen’s algorithm alone will not be 
sufficient to reduce lattice bases in higher dimensions. We may need to combine 
Seysen’s algorithm with other lattice basis reduction techniques to efficiently reduce 
large lattice bases. Alternately, it may be possible to use some heuristic technique to 
reduce the probability of reaching a lattice basis which is a local minimum of S'(A), 
or to “kick” a Seysen-reduced basis out of a local minimum. The following section 
suggests a few possible methods which may be employed when Seysen’s algorithm 
fails to .S'-reduce a lattice basis and stops at a local minimum of S'(A). 
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n 

Lio 

S(A) 

T* 

^10 

iio 

S(A') 

r * / 
^10 

# Steps 

(Seysen) 

# Steps 

(Lovasz) 

20 

82.37 

4.71 • 10 11 

118.22 

0. (20) 

20. 

0. 

757.37 

4424.0 

25 

131.95 

5.67 • 10 14 

192.93 

0. (20) 

25. 

0. 

1622.11 

12977.4 

30 

195.81 

9.28 • 10 17 

287.24 

0. (20) 

30. 

0. 

3275.30 

32718.5 

31 

205.86 

3.05 • 10 18 

308.06 

0. (20) 

31. 

0. 

3690.51 

38085.2 

32 

222.41 

1.04 • 10 19 

323.66 

0. ( 11 ) 

1695.8 

0. 

3626.00 

44770.0 

33 

237.86 

5.07 • 10 19 

348.19 

0. (5) 

86999.0 

0. 

3408.27 

53407.4 

34 

243.86 

4.48 • 10 19 

358.86 

0. ( 1 ) 

4.26 • 10 6 

0. 

3075.63 

60299.7 

35 

269.25 

6.24 • 10 2 ° 

390.02 

110.39 

1.00 • 10 8 

112.57 

3036.95 

70047.2 

36 

287.16 

5.26 • 10 21 

423.62 

138.11 

2.48 • 10 9 

141.13 

2738.39 

80836.3 

37 

291.86 

5.23 • 10 21 

438.94 

162.17 

3.20 • lO 10 

167.33 

2400.10 

91491.5 

38 

310.35 

3.17 -10 22 

469.18 

175.52 

1.07- 10 11 

182.50 

2520.02 

103624.0 

39 

329.32 

1.78 • 10 23 

500.89 

196.26 

7.45 • 10 11 

203.87 

2555.72 

118338.0 

40 

343.84 

1.82 • 10 23 

507.02 

216.88 

5.21 • 10 12 

227.10 

2240.23 

132901.0 

45 

434.01 

2.02 • 10 26 

656.67 

304.56 

4.20 • 10 15 

323.53 

2369.90 

- 

50 

541.20 

4.53 • 10 29 

833.54 

402.80 

2.33 • 10 18 

429.15 

2698.11 

_ 










































































































































2.3. WHEN SEYSEN’S ALGORITHM FAILS 


39 


2.3 When Seysen’s Algorithm Fails 

We have seen above that for random, dense lattices L with det(T) = 1, Seysen’s 
algorithm starts to break down for 30 < n < 35. As n increases beyond 35, the 
number of local minima for the 5(A) function apparently increases, and thus the 
chance that Seysen’s algorithm will reduce lattice L to I n decreases. As the number 
of local minima increases, it is increasingly likely that in the process of reducing 
lattice L the 5(A) function (where A is the quadratic form associated with L) will 
encounter one of these local minima. As described above, Seysen’s algorithm cannot 
tell whether it has reached a local or a global minimum. Thus, it stops as soon as all 
possible JA transformations cause 5(A) to increase. 

For many types of lattices, such as the sparse lattices generated by subset sum 
problems (see Chapter 3), Seysen’s algorithm has performed sufficient work by the 
time it encounters a local minimum that it is acceptable for it to stop. However, for 
many lattice reduction problems Seysen’s algorithm stops too soon. We would like 
the algorithm to be able to detect local minima and overcome them. If one considers 
the surface described by 5(A) values, local minima are “wells” or “depressions” in 
the surface which are large enough to contain all points reachable by performing one 
row move on the lattice. In this section we discuss possible techniques for “kicking” 
the reduced lattice out of these wells; methods of enhancing Seysen’s algorithm so 
that it may consider other options when it encounters a local minimum. 

There are many methods which could conceivably be applied to a lattice to move 
it out of a local minimum; we consider only some of these options. Section 2.3.1 
considers an obvious possibility, which is to consider row moves involving 3 or 4 
vectors at once (general n-moves are discussed in Section 2.4.1 below). In Section 2.3.2 
we investigate simulated annealing and rapid quenching approaches to the problem. 
Finally, Section 2.3.3 discusses using Hadamard matrices to permute the entire lattice 
basis. 
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2.3.1 Row Moves Involving Three or Four Vectors 

As initially described in [38], Seysen’s basis reduction algorithm only considers row 
moves involving two basis vectors. That is, we only consider moves of the form: 

bj 4- bj + Abi. (2.16) 

(These moves, of course, correspond to the set of transformation matrices TA.) How¬ 
ever, there is no real theoretical reason to restrict ourselves to row moves of the form 
given in Equation 2.16. We consider here possibilities for row moves of the form 

bj <- bj + Ab it + /cbj 2 , 
bj <— bj + Ab ix + /cb i2 + fi b ia , 

and their dual basis counterparts 

bj - bj + Abr, + Kb,*,, 

bj bj + Ab? + /cb,* 2 + ^b,* 3 . 

Before we begin, it should be noted that there are practical reasons for not consid¬ 
ering row moves involving more than two vectors at any given time. First, if we are 
using a greedy selection method to choose the vectors upon which to operate, more 
work is required to choose the correct n-tuple. (If we use a lazy implementation this 
is less of a concern). Second, 2-moves 3 exhibit a symmetry between the prime and 
dual lattice which is lost when we consider n-moves with n > 2. When we perform a 
TA transformation, bj <— bj + Abi and bj* <— bj* — Abj. Thus we need not explicitly 
consider operations on the dual lattice, since every dual lattice transformation has an 
equivalent transformation in the prime lattice (with i,j swapped and A multiplied by 
-1). For n-moves with n > 2, however, this duality is lost. The 3-move 

bj «— bj + Abq + «bi 2 , 

3 We use “A;-move” to designate a row operation in which multiples of k — 1 vectors are added 
simultaneously to another vector. 
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for example, has the following effects in the dual basis: 


b t 


b? 


*2 


^-AbJ, 


b? a - Kbj*. 


Thus, even if we use a lazy approach to choose which vectors to use in a row operation, 
there are many more candidate operations which must be considered since there is 
no longer any overlap between moves in the primal and dual lattices. 

Calculating the best possible values of A, k, /z, ... for a given set of basis vectors is 
also more complicated as the number of vectors involved in the move increases. As an 
example, let us consider the computation required to choose A and k for the 3-move 


bj bj + Abq + «bi 2 , 
K <- K - Ab?, 


b r 


b? - Kb?. 

12 J 


We assume without loss of generality that A < i 2 < j. Then we may represent this 
transformation by a transformation matrix, where: 

= In + A Ui u j + rdJ t2 j . 


Similar to Section 2.1.4 above, let us define 

A(n,i 2 ,i,A,K) = S((T*jy A T*j) - S(A) 

( “I” 2a; 2 ^ & T 2u^ 2 it2 &jj k 

+ 2a,- li ,- 1 a fiiJ - 4- 2A + 2a^,i 2 a jj K ^ T 

Now, we compute the partial derivatives of A(«i, i 2 ,j, A,/c) with respect to k and A, 

d 

A(?i, I 21 J 1 A, k) — 2a, 2){2 a l2 j -)- 2ctj' 2 j0tjj -(- 4<7,, 2 , 2 ctjj k T 2a,', )t2 Qjj A , 

Q 

^A(*i, * 2 tji A, k) = 2a tljt jOjj j + 2cii lt jCij'j + 2aj 1) ; 2 ajj/v + 4aj lt »jOjjA, 
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set them equal to 0, and solve for k and A: 


dA dA 


A = 


■4ai 2 ,i 2 ( a- - 2a, 1)8 - 2 a ? j ( 2a,- 2) ,- 2 a 4 - 2 j + 2a, 2J a Jl jj 

4( a *'l,»2> 


2 ) 2 ( a jj ) 2 “ 16 °ii,<i a < a ,<a( a 


* \2 
JJ/ 


2a*!,* 2 ( 2a* li j 1 a ii j + 2a * a jj 40t 1 ,* 1 Ojj ^ 2a t - 2i { 2 o, 2 j + 2o* 2 ,j a j t jj 

4( a n.*2) 2 ( a jj) 2 — l® a n»n a «2,*2( a j,j) 2 

If we wish to implement a version of Seysen’s algorithm which allows 3-moves in 
general, we must calculate six sets of (A, k) values; one set for each of 


bj 

<— bj + Abq + /cb; 2 

bj 

«- bj + Ab^ + Kb? 2 , 

(2.16) 

b ii 

< — bij -T Abj + fcbi 2 

K 

b? + Abj + /cbj* 2 , 

(2.17) 

b; 2 

<— b; 2 + Ab; x 4- fcbj 

K 2 

b? 2 + Ab? + /cbj. 

(2.18) 


Clearly including the six possible 3-moves and the eight possible 4-moves in Sey¬ 
sen’s algorithm is computationally expensive. However, there are reasons for wishing 
to do so. When Seysen’s algorithm reaches a local minimum of S(A) for some lattice 
L it is reducing, it has reached a point where any single row move increases S(A). By 
allowing the algorithm to look at 3-moves when it has run out of 2-moves to perform, 
we increase the number of configurations Seysen’s algorithm must investigate before 
it gives up. It is quite possible that one or more configurations which are 3-move 
attainable but not 2-move attainable will have Seysen measure smaller than 5(A). 
These reduced lattices would then permit the algorithm to move out of the local 
minimum and continue its reduction steps. 

We added routines to our implementation of Seysen’s basis reductions algorithm 
to implement three- and four-vector row operations. In all cases one vector was 
designated the “target” vector and integer multiples of the other two or three vectors 
were added to the target. We did not notice any significant improvement in the 
performance of the algorithm on random integer lattices of determinant 1 when these 
additional moves were allowed to occur on any iteration of the algorithm. In some 
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cases, the greedy selection scheme actually performed worse when allowed to use 
3-moves and 4-moves; usually this decrease in performance occurred because the 
algorithm “jumped ahead” too quickly by using the 3-move and would have done 
better by using a sequence of 2-moves. Three- and four-vector operations were helpful, 
however, when a lattice reduction reached a local minimum. In many of these cases 
a few 3-moves (or 4-moves, if considered) existed which would carry the lattice out 
of the local minimum and allow the algorithm to resume processing. Given these 
experiences and the increased complexity required to operation on more than two 
vectors at a time, we would suggest using n-moves when n > 2 only to move the 
lattice being reduced out of a local minimum. 

2.3.2 Simulated Annealing and Rapid Quenching 

Many combinatorial optimization problems have been successfully attacked using sim¬ 
ulated annealing, which was initially developed independently by [10, 21]. Simulated 
annealing approaches resemble local optimum algorithms, except that a random com¬ 
ponent is introduced which allows occasional “uphill” moves (moves which worsen the 
current solution to the problem according to a cost schedule). As simulated anneal¬ 
ing methods have been successfully applied to a wide variety of problems, it seems 
reasonable to consider adding simulated annealing techniques to Seysen’s algorithm 
in the hope of reducing the number of lcoal minima which cause the algorithm to stop 
before reaching a global minimum. 

Modifying Seysen’s algorithm to work along the lines of a simulated annealing 
approach would not be difficult. In the implementation of the algorithm, we simply 
need to accept row moves which increase the Seysen measure of the lattice basis. 
The probability of accepting a move which increases S(A) will depend upon the tem¬ 
perature of the reduced lattice, which starts high and decreases according to some 
cooling schedule and the reduction proceeds. It thus remains to specify the initial 
temperature of a lattice basis, the probability (as a function of temperature) of ac- 
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cepting a row move which increases .S'(A), and a cooling schedule which describes how 
temperature decreases with time/reduction steps. 

Another technique, based on physical systems, for solving combinatorial optimiza¬ 
tion problems is the rapid quenching approach. Simulated annealing slowly reduces 
the temperature of the solution, thus gradually reducing the probability of accepting 
a move to a higher energy/cost state. Rapid quenching, on the other hand, quickly 
reduces the temperature in the model, bringing the system to a minimum quickly. 
The system is then reheated to a temperature lower than the initial temperature, and 
the process is repeated. Seysen’s algorithm itself can be viewed as one iteration of a 
rapid quenching process. The heated system is the initial lattice basis, and the algo¬ 
rithm itself, by greedily reducing the Seysen measure of the lattice basis, decreases 
the temperature of the system. 

We modified our implementation of Seysen’s algorithm to simulate multiple rapid 
quenching iterations. When a lattice basis reached a minimum of the Seysen measure 
function and no single two-vector row move could decrease S(A), a randomization 
function was performed on the lattice to “heat” it and Seysen’s algorithm was subse¬ 
quently applied to the heated lattice basis. Our randomizing function chose a linear 
number of pairs of vectors (bi, bj), i ^ j, and (with equal probability) either added bj 
to or subtracted b; from bj. (This is the same randomizing operation used previous¬ 
ly to generate random integer lattices, except that we perform 0(n ) randomizations 
here instead of 0(n 2 ) as was done before.) Multiple iterations of the heating/Seysen- 
reducing process did successfully reduce lattice bases more than Seysen-reduction 
alone, although it is unclear as to how much benefit can be gained from repeated 
applications of this process. 

2.3.3 Using Hadamard Matrices to Permute Lattice Bases 

Our third and last suggestion for moving lattice bases out of local minima was sug¬ 
gested by Matthijs Coster (personal communication). Instead of randomly permuting 
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the lattice basis as random quenching does, Coster suggests using a Hadamard matrix 
H to transform a Seysen-reduced lattice basis B into lattice basis B' = B H. Recall 
that matrix H is a Hadamard matrix if it satisfies the following two properties: 

1. Each entry /i t j of H is either +1 or —1. 

2. If hi,..., h n are the rows of H, then for all i,j with i ^ j, (hi, hj) = 0. 

It is not difficult to show that if H is an n x n Hadamard matrix, then n = 2 or 
n = 0 mod 4 ([1], 2.21 Exercise 10). 

Now, consider the lattice basis B' obtained by multiplying B by a Hadamard 
matrix H. (If n ^ 0 mod 4 we may consider B and the corresponding lattice L to be 
sublattices of an n'-dimension space with n' > n,n' = 0 mod 4.) Each basis vector in 
B' is a linear combination of all the basis vectors in B , but no two B' vectors have 
similar constructions, since h; and hj differ in |n coordinates if i ^ j. The basis 
vectors in B' will have lengths « y/n times the lengths of the basis vectors in B, so 
while we obtain a good randomization of the lattice, the lengths of the basis vectors 
are still manageable. 

We should point out that the matrix H is not a linear transformation matrix; 
det(H) ^ 1 . This means that the lattice generated by B is not the same lattice 
generated by B'. However, all n-dimensional Hadamard matrices H n satisfy: 

H n H t n = n 

Thus, 

B H n H t n = nB 

and the net result of the operation is to scale B (and the associated lattice L) by a 
factor of n. Thus we can divide out the factor of n and we are left with the lattice L 
we started with. 

So, our plan of attack should be as follows. When Seysen’s algorithm stops and 
reports that lattice basis B is a local minimum of 5(^4), create B' = B H where H 
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is a Hadamard matrix. Now, Seysen-reduce lattice basis B' until a local minimum 
is reached. Then compute B" = A B'H *. Finally, Seysen-reduce basis B ", producing 
basis B"'. Bases B and B"' describe the same lattice L. Note that there is no 
guarantee that S(A f ") < 5(A), where A, A'" and the quadratic forms associated 
with B,B"' respectively. Further study is required before we may conclude whether 
Hadamard permutations provide a reasonable method for “kicking” lattice basis B. 


2.4 Extending Seysen’s Algorithm 

The description Seysen gave in [38] of his algorithm was only an outline of a lattice 
basis reduction technique. We have tried in this chapter to give both theoretical and 
empirical reasons for the choices made in implementing Seysen’s algorithm. However, 
we have only touched upon a few of the many possible combinations of techniques. 
As the next chapter shows, these choices are effective as reducing lattice bases derived 
from subset sum problems. For other lattices, their effectiveness may be in question. 
We briefly mention here some of the other possible choices for the various components 
of Seysen’s algorithm. 


2.4.1 General n- vector Row Operations 

Section 2.3.1 above discussed the possibility of extending Seysen’s algorithm to con¬ 
sider row operations involving three and four vectors at once. It is possible to extend 
these operations to encompass arbitrary k-moves where integer multiples of k — 1 
basis vectors are added to another basis vector. For fixed k, let bj be the target basis 
vector (the vector to be modified) and let bq,..., bi k _ 1 be the basis vectors to be 
added to bi in multiples of ..., A*,_i respectively. Then, after the row move, we 
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will have: 

k -1 

bj <— bj + ^mbj m , 
m=1 

b,* m <— b,* m - A m bj , for 1 < m < k — 1 . 


Now, we may solve for the new values of the quadratic forms A and A* after the 

row move has occurred. In A, the only value that changes is ajj. Its new value is: 

/k -1 \ k -1 k -1 

a j,j < a j,j d~ ( a irn,itn^m j d" 2A m (Zj m) j + 2 X m ^pO-i m ,i p • (2.20) 

\m=l / m=l m,p=l 

m^p 

In the dual lattice quadratic form A *, the values change for 1 < m < A: — 1. 

Their new values are: 


Imdm 


“ 2X ^ a Lj + for 1 < TO < k - 1 . 


( 2 . 21 ) 


If we compute A, the change in <S'(A) resulting from this row move, we find that: 


A = 


^ Ik -1 \ k -1 k -1 ^ 

r. a im,im^m j d" 2A m a i m ,j d~ 2A m ^pdi m ,i : 

V.m=l / m=l m,p=l 

m^p 


fc-1 


d' C Wm 2A m a im j) , 

m=1 

A :—1 fc -1 

^ = 2 ^ m a j,j a im i im + 2 ^ A m — a Wm a t m ,j) 


m=l 
A—1 


m=l 


d” 2 ^m^p a i m ,i p a jj- 

m } p~l 

rrr£p 

Thus, for all 1 < m < k — 1 we have: 

d 


k -1 


( 2 . 22 ) 


(2.23) 


3A, 


^ 4A m a,- rot i m aj i j + 2(ajja 8 - m j A 2 ^p a i m ,i P a i,i- (2.24) 


p=i 

p^m 


We may thus compute formulae for all A m for 1 < m < k — lifwe solve the si¬ 
multaneous system of equations obtained by setting the k — 1 derivatives defined in 
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Equation 2.24 equal to zero. Of course, we still must solve the problem of finding 
optimal integer values for the A m . For the 2-move case we showed that |"A r J was the 
best integer choice for A. It is not clear that rounding will work in the k- move case. 
The real values derived for the A m may only indicate some range of integer values 
which need to be searched. 

2.4.2 Alternate Selection Criteria 

Seysen’s original implementation of his basis reduction algorithm used a “lazy” ap¬ 
proach for choosing pairs of basis vectors to reduce. Pairs of integers 1 < i,j < n 

were searched in lexicographic order until a suitable row move involving bj and bj 
was found. We have presented above empirical evidence which favors a “greedy” ap¬ 
proach, even when the extra computation time required to implement the “greedy” 
method is considered. 

Selection methods other than the “greedy” and “lazy” approaches were not con¬ 
sidered in our experiments, but are certainly possible. For example, in addition to 
taking into account the reduction in S'(A) which will result from a row move, we 
might also wish to consider the other row moves which will be blocked by performing 
this move. That is, if Abj is added to bp the potential 5(A) reductions of all other 
row moves which involve either bj or bj will be modified. Perhaps we should choose 
row moves so that the moves they block have minimum 5(A) reduction potential. 
We could combine this idea with the “greedy” method; selecting a row move with 
the greatest difference between the amount it reduces S(A) and the average 5(A) 
reduction of all the moves it blocks. 

2.4.3 Alternate Choices of A 

In Section 2.2.2 above we looked at the effect of placing restrictions on the possible 
values of A on the performance of Seysen’s algorithm. In particular, A was allowed 
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either to be any integer value, or to only take on the values il. We found that the 
algorithm worked slightly better if row moves were chosen with A 6 Z but only re¬ 
moves with A = ±1 were actually performed, probably because the changes made 
to the lattice basis on any single row move are smaller. We pay for the improved 
performance, however, though an increase in the running time of the overall algorithm, 
as it takes more row operations with the restriction in place to add a large integer 
multiple of one basis vector to another basis vector. 

As an example of other possible methods of choosing or restricting A values, con¬ 
sider the following set of restrictions: 


• When choosing a pair of vectors for a row move, A may take on any integer 
value. 


• When the row move is actually performed, if A > 0 use the largest power of two 
strictly less than A (unless A = 1, in which case A should be used). If A < 0 
use the smallest power of two strictly greater than A (again, unless A = — 1, in 
which case —1 should be used). 


We may abbreviate this set of conditions as (Z,2 fc ). What we are doing is computing 
the best possible value of A, but instead of performing one row move to compute 
bj <— bj + Abi, we perform a logarithmic number of moves. In this way we might 
be able to combine the benefits of the (Z,Z) approach (fast running time) and the 
(Z, ±1) approach (better overall performance). This approach has not been tested, 
and judging from the relative differences noticed between the (Z, Z) and (Z, ±1) cases 
is not likely to produce very large changes in reduction performance. However, it is 
an example of other possible A restrictions which could be tried. 
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2.4.4 Alternate S(A ) Functions 

Section 2.1.3 above gave reasons for using functions of the product terms ||bj|| |jbj ||. 
In particular, the function 

5(A) = f:«M<i = El|bi|| 2 ||br|| 2 
2 = 1 2 = 1 

was selected as the Seysen measure function because it yielded a closed form solution 
for the optimal value of A given i and j. However, other functions could certainly 
be employed as the method of comparing different lattice bases. In this section we 
briefly describe how Seysen’s algorithm would have to be modified to accommodate 
another measure function. 

We restrict our attention to versions of Seysen’s algorithm which use only row 
moves involving two basis vectors (i.e. bj <— bj + Ab;). Recall that the formula in 
Equation 2.9 for the optimal choice of A was derived by maximizing the change in the 
Seysen measure function caused by a row move involving two particular basis vectors. 
In the Seysen measure function is changed, the only direct impact it will have upon 
the operation of the algorithm is that the optimal value of A for basis vectors (bj, bj) 
will be computed in a different manner. 

In [38] Seysen mentions two possible replacements for the S'(A) function: 

Si(A) = f[«2,2< i = f[ UbilHlKH 2 , 

i=i i —i 

S*{A) - E \/«m<, 2 = E llbillllbi ||. 

2 = 1 2 = 1 

Replacing 5(4) with 5 1 (A) implies that our choice of A must minimize: 


I . \ 


n 




. m=1 


( a 2,2 a 2,2 a j,j a j,j a i,i a i,i a j,j a j,j) ' 


) 


Solving = 0 yields an extremely complex expression for A. Similar results occur 



2.4. EXTENDING SEYSEN’S ALGORITHM 


51 


when we try to substitute £ 2 (.<4) for 5(A). In both cases, no simple closed-form 
solution for A exists as was the case with 5(A). 

It may be still be possible to utilize measure functions in Seysen’s algorithm with 
no simple closed-form solution for A if we are willing to sacrifice some performance. 
If the range of possible integer A values is bounded, for given i and j we can compute 
A(i,j, A) for all possible A values in the permitted range. The A which provides the 
greatest change may then be selected. The cost of this procedure is that we can 
no longer guarantee that the maximal A for a pair of basis vectors (bi,bj) may be 
found in constant time. If the range of A values to consider is usually small, then we 
will probably notice little more than a linear slowdown in the running time of the 
algorithm. For large ranges of possible A values, further heuristics might be applied, 
such as only considering A values which are near a power of two. 

In our experiments we noticed that large A values tended to occur only during 
the first few row reduction performed by Seysen’s algorithm. After this initial burst 
in reduction of 5(A) row moves tended only to involve small integer A values; it was 
quite rare to find |A| > 10. If similar conditions occur for the lattice bases in question, 
it is probably reasonable to use a more complex measure function than 5(A) and use a 
small exhaustive search over a bounded range of possible A values to find the optimal 
A coefficient for a row move. 






Chapter 3 


Solving Subset Sum Problems 

3.1 Introduction 

In the previous chapter we discussed the theoretical and empirical implications of Sey- 
sen’s basis reduction algorithm. As Chapter 1 pointed out, many problems in discrete 
mathematics may be reduced to problems involving lattices and basis reduction. In 
this chapter we use Seysen’s algorithm to solve subset sum problems. The subset sum 
problem in the general case is NP-complete, and many knapsack-type cryptosystems 
have been suggested which depend on the difficulty of solving subset sum problems. 
We show below that Seysen’s algorithm, when used in conjunction with the LLL al¬ 
gorithm and other techniques, allows us to solve a large class of subset sum problems 
in polynomial time. 

We begin our analysis with the definition of a subset sum problem. 

Definition. Let A — {ai,...,a n } be a set of positive integers (the weights). Let 
A' C A be some subset of A, and let s be the sum of the elements of A'. Then 

n 

s = ’Yf l t i a i , for ti € {0,1}, 1 < < n. (3.1) 

*=i 

The subset sum or knapsack problem is to find, given the set A of weights and the 
sum s, some subset A' of A. Equivalently, one may find a set of values for the 0-1 
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variables ei,..., e n , where e» = 1 if and only if a{ E A'. 

The subset sum problem is known to be NP-complete [15]. In the INSTANCE- 
QUESTION format often used to phrase NP-complete decision problems, the subset 
sum decision problem would be described as follows: 

INSTANCE A set of positive integers A = {ai,..., a n } and an integer s. 

QUESTION Does there exist a subset A' of A such that the sum of the 
elements of A' is s? 

Clearly, if we can solve subset sum problems in polynomial time, we can answer 
the subset sum decision problem in polynomial time. The converse is also true; we 
can ask an oracle which answers the subset sum decision problem in polynomial time 
whether there is a solution to the subset sum problem with weights a 2 ,... ,a n and 
sum s — a\. If there is such a solution, then there exists a solution to the original 
problem that included a\ (i.e. we can set e\ — 1). If the oracle say no such solution 
exists, then a\ cannot be in the subset which sums to s, and we know that e\ = 0. 
We can then recurse and determine e 2 , e%, ..., e n in sequence. 

Many public-key cryptosystems have been proposed with the difficulty of solving 
subset sum problems as the basis for their security. (See [7, 8, 13, 31] for surveys of 
this field.) Almost all of these cryptosystems have been shown to be insecure; the 
Chor-Rivest one [11] is perhaps the most widely known system which has not yet been 
broken. The majority of the attacks on knapsack-based cryptosystems have involved 
discovering the secret information hidden in the weights which allows the receiver A to 
decrypt the message quickly. However, there have been two independent attacks, one 
due to Brickell [6] and one due to Lagarias and Odlyzko [26], which attempt to solve 
all subset sum problems of a certain type, independent of the method in which the 
weights were chosen. These methods (and the newer result in [12]) depend in theory 
only on the density of the subset sum problem to be solved. In practice, however, the 
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success rate of these methods is bounded by the performance of the basis reduction 
technique used in the attack. 

Section 3.2 below outlines currently known methods for solving subset sum prob¬ 
lems, and describes a new method which significantly increases the class of problems 
which may be attacked in practice. Section 3.3 discusses current methods for actually 
solving a specific subset sum problem, and the limits of these methods. In Section 3.4 
we show how to use Seysen’s algorithm in conjunction with multiple versions of the 
LLL algorithm and other heuristics to solve a larger class of subset sum problems. 
Section 3.5 presents empirical results obtained by solving a large number of subset 
sum problems using Seysen’s algorithm. These results give experimental evidence 
that it is possible to solve a much larger class of subset sum problems in polynomial 
time than was previously thought possible. 

3.2 Theoretical Bounds on Solving Subset Sum 
Problems 

The majority of attacks on knapsack-based cryptosystems exploit the specific con¬ 
structions of the cryptosystems. Two algorithms have been proposed, however, which 
depend only on the properties of the subset sum problem and not on any specific 
method of construction. These algorithms, one by Brickell [6] and one by Lagarias 
and Odlyzko [26], show that almost all low-density subset sum problem may be solved 
in polynomial time. The density d of a set of weights <q,..., a n is defined by 

d=-. ---• (3.2) 

logo max a* 

For d > 1 there will in general be many subsets of weight with the same sum s, so 
from an encryption-decryption point of view we are interested in subset sum problem 
instances with d < 1. The Brickell and Lagarias-Odlyzko algorithms show that it is 
possible to solve almost all subset sum problems with d sufficiently small. 
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The Brickell and Lagarias-Odlyzko attacks reduce the subset sum problem to 
the problem of finding the Euclidean-norm shortest nonzero vector in a lattice. As 
was mentioned in Section 1.1 above, finding short vectors in lattices may be very 
hard in general. The theoretical worst-case bounds for the LLL algorithm and its 
variants are not encouraging, and no bound currently exists for Seysen’s algorithm. 
However, these techniques tend to perform much better in practice than in theory; the 
performance of Seysen’s algorithm on the Korkin-Zolotarev test lattice with 9 = 0.4 
(Section 2.2.3) is one such empirical example. Thus, it seems important to separate 
the efficiency of lattice reduction and finding short nonzero vectors from the difficulty 
in reducing subset sum problems to lattice reduction questions. 

We consider a Euclidean-norm lattice oracle (or lattice oracle for short) that, when 
given a lattice basis as its input, with high probability finds in polynomial time the 
Euclidean-norm shortest nonzero vector in the given lattice. We do not know how 
to construct such an oracle, but it might be possible to do so. Data provided in 
[26, 33] show that at low densities the LLL algorithm behaves as a lattice oracle, and 
our results in Section 3.5 below show that this is also the case for a combination of 
the Seysen and LLL algorithms, even for significantly larger and denser subset sum 
problems. Given the existence of a lattice oracle, the analysis in [26] shows that 
it is possible to solve almost all subset sum problems of density d < 0.6463... in 
polynomial time. Recently, Coster, LaMacchia, Odlyzko and Schnorr [12] and Joux 
and Stern [19] independently demonstrated via different techniques that this bound 
could be improved to d < 0.9408. In fact, if we assume the existence of a sup-norm 
lattice oracle instead of a Euclidean-norm lattice oracle, [12] showed that the density 
bound then becomes d < 1. 

The Lagarias-Odlyzko attack proceeds as follows. Let {cq,..., a n } be a set of 
weights with 0 < a,- < A for some positive integer A and for all 1 < i < n. Let 
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e = (e x ,..., e n ) G {0,1}", e ^ (0,0,..., 0) be fixed and depend only on n. Then 

n 

S = ^2 e i a h f° r e * £ {0, 1}, 

i =1 

is the sum of the subset of weights, where a,- is in the subset if and only if e; = 1. 
Now, define basis vectors b x ,..., b n+x as follows: 

b x = (l,0,...,0,JVai), 
b 2 = (0,1,..., 0, IVa 2 ), 


b n = (0,0, ...,l,jVa n ), 


b„+i = (0,0, ...,0,iVs), 

where TV is a positive integer > y/n. Let L be the lattice defined by the basis vectors 
bi,...,b n+1 . That is, 


L = | Zibi: Z{ G Z, for 1 < i < n + 1 j . 

Notice that lattice L contains the vector e = (e x , e 2 ,.. •, e n , 0), the solution vector to 
the subset sum problem, since 

n 

e = ^2 e t b; - b n+x . 

i =1 

Let P denote the probability that there exists another vector x G L such that 
||x|| < ||e|| and x ^ {0,e, —e}. The simplified analysis of the Lagarias-Odlyzko 
attack presented in [14] shows that this probability is bounded: 

P < n (2nJ\k+ l) for c 0 = 1.54724... (3.3) 


Thus, if the bound on the size of the weights A = 2 cn with c > Co, lim P — 0. If the 

n—+oo 

density of a subset sum problem is less than 0.6463..., then 


log 2 max a, 

l<i<n 


< 0.6463.. 


__ _ __ ryn/0.6463., 

max a; > Z 1 

l<i<n 


A > 2 C °". 
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Thus, all subset sum problems with density < 0.6463 ... could be solved in polynomial 
time, given the existence of a lattice oracle. 

Recently, two independent improvements [12, 19] to the Lagarias-Odlyzko attack 
have been developed, both of which increase the density bound to d < 0.9408.... 
The modification suggested in [12] is to replace the b n+ x basis vector in L with 

K + i = (bb---,lNs). 

Let L' be the lattice spanned by the vectors bi,..., b n , b„ +1 . Lattice L' does not 
contain the solution vector e, but it contains a similar vector e': 

e' = (ei,..., e' n , 0), where e- = e t - - 

We know that e[ G { —1,|} for 1 < i < n since e t - G {0,1} for 1 < i < n. Then 
||e'|[ 2 < \n independent of the number of e;’s which are equal to 1. 

Using lattice L' we are now interested in the probability P' that there exists a 
vector x' G L' such that: 

Ill'll < lie'll < \y/n, 

\oA) 

x' ^ {0,e',-e 7 }. 

Utilizing similar techniques to those in [14, 26, 28], [12] shows that the probability P' 
is bounded above by: 

P' < n (iny/n + 1^ for c' 0 = 1.0628- (3.5) 

This bound is similar to that in Equation 3.3 above. Since l/c' 0 = 0.9408..., any 
subset sum problem with density d < 0.9408... may be solved in polynomial time, 
given the existence of a lattice oracle. 

3.3 Previous Empirical Methods 

Along with their theoretical results, Lagarias and Odlyzko [26] presented in 1985 
the results of the first empirical attacks on general subset sum problems. Their 
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method was to apply a multiprecision version of the LLL algorithm to the basis L 
presented in Section 3.2 above (b n+ i = (0,0,0,... ,0, s)). Lagarias and Odlyzko set 
the LLL parameter y — 1, which they said yielded better results than y — 0.75 but 
tripled execution time in practice. Also, five random orderings of the basis vectors 
for lattice L were tried, since different initial orderings yield different LLL-reduced 
bases. Experiments were conducted on various subset sum problems with n < 50 and 
densities d between 0.5 and 0.875. Figure 3-1 graphically shows the performance of 
the Lagarias-Odlyzko method. For each tested value of n, the labeled curve connects 
(density,success rate) points in the plot obtained by attempting to solve subset sum 
problems of size n. 

For n < 26, the LLL algorithm appeared to function almost as well as a square- 
norm lattice oracle; all subset sum problems with density d < 0.6408... and n < 26 
were solved when LLL was used on five random permutations of the basis vectors. 
However, performance degrades quickly as n grows above 30. For n = 40, Lagarias 
and Odlyzko were able to solve all attempted subset sum problems only for density 
0.5. At n = 50, only two-thirds of the attempted density 0.5 problems were solved. 

In [33] Radziszowski and Kreher reported the results of extensive attempts to solve 
subset sum problems. Their reduction algorithm, based on LLL, differed from that 
used in Lagarias-Odlyzko in two important ways. First, Radziszowski and Kreher 
modified the LLL algorithm to reduce the number of required multiprecision calcula¬ 
tions. In essence, instead of running LLL once on a lattice with numbers k bits long, 
they ran LLL m times with numbers k/m bits long (a divide-and-conquer approach). 
This modification sped up their algorithm and allowed them to attack larger subset 
sum problems than Lagarias and Odlyzko (larger here means greater value of n). 

The second change Radziszowski and Kreher made to the Lagarias-Odlyzko attack 
was to use another algorithm in conjunction with LLL to find short vectors in the 
lattice. This second algorithm, called Weight-Reduction, searches for pairs of vectors 
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(bi,bj) in the lattice for which: 

||bi + ebj|| < max{||bi||, ||bj||}, for e = ±1. (3-6) 

When Weight-Reduction finds such a pair, it replaces the larger of bj and bj (in terms 
of square-norm) by the sum bi + ebj. Equation 3.6 is satisfied by a pair («, j),i ^ j 
if and only if 


max{||b i || 2 ,||b j || 2 }<2-||(b i ,b j )||. 

The Weight-Reduction algorithm can easily be implemented in 0{n 2 ) time to search 
for all pairs of vectors (bpbj) which satisfy Equation 3.7. 



Figure 3-2: Radziszowski-Kreher Inner Loop 

Radziszowski and Kreher alternated calls to the modified LLL algorithm with calls 
to Weight-Reduction and a sorting procedure to improve the search for short vectors 
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in the lattice. Figure 3-2 is a flowchart-like representation of the process. The initial 
basis L is LLL-reduced (y = 0.99 for all invocations of LLL). The output of LLL is then 
Weight-Reduced, sorted by length, and then fed back into LLL. This process continues 
until some termination condition is met. In theory, “no decrease in ||bi||” could 
be used as the termination condition. In their experiments, Radziszowski and Kreher 
set an explicit bound on the maximum number of loops which could be performed. 
The output of this iterative process is then Weight-Reduced one more time. 

One important feature of the Radziszowski-Kreher loop is the insertion of the 
sorting phase before LLL is run each time. By sorting the vectors in the lattice 
by length, the length of the shortest vector in L is guaranteed not to increase by 
application of LLL. If L is a lattice output by the sort procedure, then vector bi in 
L is the shortest vector in the lattice. Now, LLL can replace bi with some vector 
bj only if HbjH < ||bx||. (This does not hold for bi in general but is true for b x .) 
Thus, we are guaranteed that the shortest vector in L before LLL is applied will not 
disappear from the basis unless an even shorter vector is found. 

Figure 3-3 shows graphically the results of the experiments carried out in [33]. 
For these experiments, only 15 loop iterations were allowed (9 if all the vectors in the 
basis were of length < \n). For values of n < 34, this algorithm was able to solve all 
attempted subset sum problems of density d < 0.654, which is above the Lagarias- 
Odlyzko bound. As n increases from 42 to 98 the density at which all attempted 
problems were solved decreases. As was reported in [26], the effectiveness of LLL as 
a lattice oracle drops off as the size of the lattice grows. Radziszowski and Kreher’s 
results show however that it is possible to increase the effective range of LLL by 
combining it with other heuristic techniques. 

The Radziszowski-Kreher algorithm was the best known method to date for solving 
subset sum problems. In the following section, we show how to to combine the 
Seysen basis reduction algorithm with the LLL algorithm, other heuristics, and the 
theoretical improvements in [12] (Section 3.2) to greatly extend the range of subset 
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sum problems which may be empirically attacked and successfully solved. 


3.4 Using Seysen’s Algorithm to Solve Subset Sum 
Problems 

In both [26] and [33] the LLL algorithm was used almost exclusively to perform 
the required lattice basis reductions. Given the existence of Seysen’s basis reduction 
algorithm and knowledge as to how it performs relative to LLL, it is natural to wonder 
how Seysen’s algorithm could be applied to solving subset sum problems. We know 
from some of the comparisons performed between LLL and Seysen’s algorithm in 
Chapter 2 that Seysen’s algorithm tends to perform fewer row operations than LLL. 
If we used Seysen’s algorithm then as part of an attack on subset sum problems, 
we should be able to perform more reduction operations (or other work) without 
increasing the overall execution time of our method. 

There are also theoretical and empirical reasons to suspect that Seysen’s algorithm 
alone will not perform well on subset sum problems. Seysen’s algorithm was originally 
suggested for simultaneously reducing a lattice and its dual lattice, not for finding the 
shortest vector in a lattice. If a lattice and its dual are of comparable size, the Seysen 
algorithm is not likely to perform row operations that generate short vectors in one 
lattice if that move increases the lengths of vectors in the dual lattice. If the vectors 
in the dual are significantly larger than those in the prime lattice, Seysen’s algorithm 
may actually produce vectors larger than those originally input. Also, the algorithm’s 
operation is not dependent on the ordering of the basis vectors, which means that 
the benefits gained by running LLL multiple times on different randomizations of the 
subset sum lattice disappear. 

These predicted benefits and deficiencies of the Seysen algorithm suggest that 
the best place to use Seysen is at the beginning of our attempt to solve a subset 
sum problem. The initial subset sum lattice L has primal basis vectors which are 
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much larger than those in the dual; this favors the generation of shorter vectors in 
L and longer vectors in L*. Also, as Seysen tends to perform many fewer row moves 
than LLL to reach the same reduced lattice basis, using Seysen’s algorithm initially 
will reduce the total number of multiprecision row moves as compared to using only 
LLL. Once Seysen’s algorithm stops at a local minimum, we can use LLL and other 
techniques to further reduce just the primal lattice L and look for short vectors. 





Seysen-Reduction 


' 



GCD-Reduction 


' 

i 


Seysen-Reduction 


* 

f 



Figure 3-4: Overview of Algorithm SL 

Figure 3-4 shows the basic outline of the SL (Seysen-Lovasz) algorithm. We 
start with the initial lattice basis suggested in [26]. If the weights of the subset sum 
problem are sufficiently large, then the lengths of the basis vectors in L will be much 
greater than the corresponding lengths of the basis vectors in the dual lattice L*. We 
apply the Seysen algorithm to the L , L* self-dual pair of lattice bases, using a greedy 
algorithm (Section 2.1.2) to choose at each step the pair of basis vectors (b;,bj) to 
be reduced. The measure function S'(A) is the one suggested by Seysen [38]: 

•s , (-4) = EI|b i |i 2 iib:ii 2 . 

i=1 

When computing the change in S(A) for each pair of vectors (bp bj), i ^ j, A is always 
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set to its maximum value of: 



Let Ls t be the lattice basis which is produced by this (first) application of the Seysen 
algorithm. 

Recall that the goal of the SL algorithm is to find a basis vector which describes 
the solution to the given subset sum problem. In [26, 33] the desired solution vector e 
was always of the form (ei,..., e„) with e 4 - equal to either 0 or a fixed constant k. We 
will call all such vectors e Type-I solution vectors, to designate them from another 
class of solution vectors which will appear below. At the completion of the first Seysen 
reduction in the SL algorithm, it is possible that lattice basis Ls 1 contains a Type-I 
solution vector. If so, then we are finished, and Algorithm SL halts. We assume that 
Ls 1 does not contain a Type-I solution vector, and continue with the next stage of 
the algorithm. 

The algorithms of [26, 33] searched for short vectors of the form (ei,... , e n , e n+ i = 
0) in the LLL-reduced bases, where e; {0,«} for 1 < i < n and /c 6 Z is any fixed 
integer. One of the problems with these methods is that often the short vectors 
produced by LLL reduction had e n +i ^ 0; that is, the sum X) e i a i was n °t °f the 
form 5 + yt and did not describe any relation involving the target subset sum. One 
method of reducing the probability that LLL (or Seysen) will include a vector of 
this form in the reduced basis is to scale the a,- and s by some constant factor N, 
which increases the length of any vector having a nonzero e n+ i by about y/N if N 
is sufficiently large. However, this approach has the drawback that it increases the 
size of numbers which are already quite large and require multi-precision arithmetic. 
We suggest another method for elimination from consideration all lattice vectors with 
e n+ i ^ 0: GCD-reduction. 

The idea behind GCD-Reduction is to perform row moves on the lattice basis so 
that the entire (n + l) th column contains exactly one nonzero element, the GCD of the 
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weights ai,... , a n . If (for example) vector b n+ i contains the GCD, we can remove 
b n+ i from the basis and remove the (n + l) st column completely from the lattice 
basis. This reduces what was an (n + l)-dimension lattice to an n-dimension lattice, 
and guarantees that any lattice vector generated by reducing this basis would have 
had its (n + l) st component equal to 0 in the (n + l)-dimension space. 

Implementing GCD-Reduction is easy. The basic algorithm we use was described 
by Brun [9]. Basis vectors are sorted in order of decreasing 6,- >n +i. Then b 2 is re¬ 
peatedly subtracted from bi until & 2 ,n+i > &i,n+i* The vectors are then resorted in 
order of decreasing 6;,„+i and the process loops. In reality, only vector bi needs to 
be inserted into the previously generated order. This can be performed very fast by 
using an auxiliary pointer array; the actual basis vectors don’t even need to be moved 
around. Eventually, the only nonzero element in the (n + l) st will be &i, n +i! at which 
point both the vector bi and all 2 < * < n + 1 can be removed from the lattice. 

We apply GCD-Reduction to the lattice basis Ls x output by the first application 
of Seysen’s algorithm. We do not use GCD-Reduction on the original lattice basis 
L because the resulting lattice basis matrix would be quite dense (which means the 
Seysen algorithm will take longer to run) and the vectors in the dual lattice would also 
be much larger. After GCD-Reduction is applied to basis Ls 1 (yielding output basis 
Lg in n dimensions) we again search for the presence of a Type-I solution vector. 
Assuming we do not find such a vector, the Seysen algorithm is applied to Lg to 
reduce the lattice to a local minimum of the measure function S'(A). Call this lattice 
Ls 2 i the output of the second Seysen application. 

For small values of n (n < 20), it is possible that Ls 2 will contain a Type-I solution 
vector. For larger values of n and sufficiently high densities d, however, it is unlikely 
that the Seysen- GCD-Reduction-Seysen portion of the algorithm will have found the 
desired vector. Thus we begin the second portion of the algorithm: the LLL phase. 

Figure 3-5 shows the first level of detail in the LLL-phase of Algorithm SL. The 
input to this phase of the algorithm is the lattice basis Ls 2 which was the result of 
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Figure 3-5: Algorithm SL: LLL-Phase 
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the second application of Seysen reduction. We now take advantage of the theoret¬ 
ical improvements in the density bound given in [12] and Section 3.2 above. (The 
introduction of the “one-half” vector is delayed until after the Seysen stage of Al¬ 
gorithm SL has completed to increase the sparsity of the lattice bases which are 
Seysen-reduced.) We extend the R n lattice described by Ls 2 into a lattice in R n+1 by 
adding an (n -f l) st component to each of the basis vectors bi,..., b n : 


old b; — ipi.lt bi, 2 t • • • t bi,n ) 


new bj = (&;,!, b it2 , ■■■, &»>+i = hi) 

j =i 


The new (n + l) st component of each basis vector is simply the sum of the first n 
components. To complete the extension, we need to add one more vector to the lattice 
basis. New basis vector b n+1 is defined as follows: 


bn+l ( 2 ’ 2’ ' ’ ' ’ 2’ 2^ 1)) 

n terms 

For b n+ i, the (n + l) st component is the sum of the first n terms minus The — | 
correction is needed because otherwise the (n + l) st column of the basis would be 
dependent. Let L m designate this extended lattice. 

With the introduction of vector b n+ i we have also introduced a second class of 
solution vectors. If e is a Type-I solution vector for a subset sum problem with 
elements in the subset, then there now exists a Type-II vector in the lattice of the 
form 


e = (ei,..., e n , e n+1 ), for e,- G {-§, |}, 1 < i < n + 1. 

Also, if e was a solution vector of the n-dimensional lattice, in the extended lattice 
T; n vector e will be transformed into e' = (ei,..., e„, e„ +1 ) where e n+ i equals k times 
the number of elements in the subset which sums to s. Thus, during the LLL phase 
of Algorithm SL, we search for both extended Type-I solution vectors and Type-II 
solution vectors. 
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The structure of the LLL phase of Algorithm SL is based upon the successful 
attacks of [26, 33]. Lagarias and Odlyzko showed that using LLL multiple times on 
random orderings of the input basis improved the chances of finding a solution vector. 
Thus, Algorithm SL initially tries to reduce lattice Li n , and upon failure randomizes 
L m and tries again. This process continues until a total of tti reduction attempts 
have been made (L- in and -K\ — 1 random orderings of L m ). In [26] 7Ti = 5, but there 
is nothing special about that particular value. 

Recall from Section 3.3 above that the Radziszowski-Kreher approach, while set¬ 
ting 7Ti = 1, repeatedly called LLL in conjunction with Weight-Reduction and Sort 
(Sort just sorts the basis vectors by length.). They showed that repeated iterations 
of the LLL- Weight-Reduction-Sort loop shown in Figure 3-2 yielded better results 
than a single application of LLL. Algorithm SL incorporates this approach, applying 
a number of iterations of LLL- Weight-Reduction-Sort before giving up and trying a 
new randomization of L ia (Figure 3-5). Each ordering of L m passes though up to 7r 2 
iterations of the loop before Algorithm SL gives up. As stated above, in [33] 7r 2 = 15, 
or 7T 2 = 9 if all the vectors in the lattice were of length < ^n. Again, there is nothing 
special about these particular values. Theoretically, one could choose 7Ti and 7r 2 to be 
quite large. Any practical algorithm, though, would probably want to use reasonable 
constants for both 7Ti and 7r 2 . 

Although both [26] and [33] run the LLL algorithm with the parameter y « 0.99, 
Lenstra, Lenstra, and Lovasz show that y may be any value in the range \ < y < 1 
[27]. For y & LLL will only exchange vectors in the lattice (Equation 1.2) for 
relatively large differences between vectors b; and bi + i. As y —*■ 1, the amount 
of improvement required to trigger the LLL exchange step decreases, and LLL will- 
swap vectors and continue running for minimal improvement. Thus, larger values 
of y will likely lead to improved results, but LLL will also take longer to run. This 
suggests that instead of running the entire LLL algorithm with y = 0.99, as was done 
previously, a version of LLL which started with y ~ | and adjusted it upwards as the 
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Figure 3-6: Algorithm SL: LLL-Loop Structure 
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algorithm ran might work just as well but require fewer row operations [18]. 



Figure 3-7: Algorithm SL: LLL(x) Internal Structure 


Algorithm SL uses a version of the LLL algorithm which varies the parameter y 
among a finite number of values. Figure 3-6 shows what happens each time SL calls 
LLL-Loop. For subset sum problems with n < 90, each call to LLL passes through 
six stages (i.e. six distinct values for y) 1 . Upon entry to each stage, the y parameter 
is updated, and LLL-reduction is performed upon the lattice basis using the new y 
value (see Figure 3-7). Let L VOt ; n and L yo)OUt represent respectively the lattice bases 
input to and output from an application of the LLL algorithm with y = yo . We now 


1 Three additional stages are used for subset sum problems with n > 90 to help reduce error. 
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compute the boolean value of the following expression: 

( min{||b;|| £ L yo<oxl t } <C min{||bi|| G -hyoiin}) 

OR( (min{||bi|| G T^out} = min{||bi|| G L y0iin }) ( 3 - 7 ) 

AND (max{||bi|| G Ty 0 ,out} < max{[|bi|j G -hy 0 ,in})) • 

(This decision point is represented by the diamond in Figure 3-7.) If the boolean 
value of Equation 3.7 is TRUE, then the Weight-Reduction procedure is applied, the 
Sort procedure is invoked, and the output is input again into LLL with y = yo- Thus, 
we perform LLL- Weight-Reduction-Sort loops until either the length of the shortest 
vector in the lattice has increased after LLL-reduction, or if the length of the shortest 
vector has remained constant and the length of the longest vector has not decreased 2 . 

By using the recursive structure in Figure 3-7 and the termination condition rep¬ 
resented by Equation 3.7, we attempt to have the LLL algorithm perform as many 
reduction steps as possible with small values of y 0 . We delay increasing the value 
of y until LLL fails to make any improvement in the length of the largest vector in 
the lattice basis. In this way we restrict the set of row swaps and moves which LLL 
will consider, which improves the running time of the algorithm. Initially, for y small 
(say y < 0.75) LLL will only consider row moves which will “significantly” reduce 
the current lattice basis. Later, for values of y « 1, LLL will consider any row move 
which reduces the lattice. 

For n < 90, six specific y values are used: 0.2578125, 0.625, 0.75, 0.875, 0.9375, 
and 0.9921875. These values were chosen for two reasons. First, all of the numbers 
have exact fractional binary representations in double-precision floating point. This 
meant that error would not be introduced into the LLL algorithm by performing 
arithmetic operations on y. Second, experiments with small subset sum problems 
(n < 24) showed that the number of row moves performed by LLL for each of the four 
middle values were approximately equal. That is, work was evenly distributed across 
2 The loop termination condition in Equation 3.7 is strictly a heuristic. Other terminating con¬ 
ditions could certainly be used. 
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all intermediate stages represented in Figure 3-6. (The endpoint values 0.2578125 
and 0.9921875 were fixed). For n > 90, three additional values were used: 0.34375, 
0.4375, and 0.53125. These values evenly divide the interval [0.2578124,0.625] into 
four equal pieces. As discussed in Section 3.5 below, our initial implementation of 
the LLL algorithm did not work correctly for lattice bases with n > 90 because of 
rounding error introduced into the computations. We were able to reduce the effects 
of rounding error to acceptable levels by modifying the LLL algorithm itself and by 
introducing the 0.34375, 0.4375 and 0.53125 stages into the algorithm. 

We have detailed the operation of Algorithm SL, a combination of the Seysen and 
the Lenstra-Lenstra-Lovasz basis reduction algorithms which also utilizes the GCD- 
Reduction, Half-Vector, Weight-Reduction, and Sort heuristics. In the next section 
we present the results of our experiments with this algorithm and show how the 
combination of these techniques greatly increases the range of subset sum problems 
which may be solved. 

3.5 Empirical Tests Using Algorithm SL 

We now present the results of our experimental attempts to solve subset sum problems 
using Algorithm SL as describes in Section 3.4 above. Following the tabulated results 
of Radziszowski and Kreher, we attempted to solve random subset sum problems of 
size n, where: 

n € {42,50,58,66,74,82,90,98,106}. 

For each value of n, a set of b values (representing the number of bits in the binary 
representation of the weights af) was chosen. Algorithm SL was then run on a number 
of randomly generated subset sum problems for each pair ( n,b ). Where possible, 
values of b were chosen to coincide with values used in [26, 33]. 

Algorithm SL was implemented in FORTRAN and utilized Bailey’s multiprecision 
floating-point package [4]. The Seysen reduction algorithm stored all lattice bases in 
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multiprecision floating point and used both multiprecision and double-precision float¬ 
ing point operations. The size of the multiprecision floating-point representation was 
changed for each value of n so that b was at most one-half the size of the representa¬ 
tion in bits. In many cases significantly more bits were available in the representation 
to help reduce rounding error. 

The LLL portion of Algorithm SL was also implemented in FORTRAN but used 
only integer and double precision operations. One advantage of the Seysen-LLL 
structure of Algorithm SL is that by the time the LLL stage is reached, the lattice 
basis to be reduced contains only relatively small integers. In all cases investigated, 
the input lattice basis to the LLL phase of the algorithm did not contain any single 
coefficient larger than 2 16 . This meant that we could safely store the basis vector 
coefficients as 32-bit integers, and we could use integer and double-precision floating 
point to carry out all calculations necessary to run LLL. Avoiding multiprecision row 
moves during the LLL phase greatly decreases the execution time of Algorithm SL. 

For values of n < 90, these implementations of the Seysen and LLL algorithms 
were sufficient for our purposes. However, for n 6 {90,98,106}, rounding error started 
to significantly alter the operation of the algorithm. Recall that the LLL algorithm 
stores the values of the coefficients as rational numbers; our implementation of 
LLL used double-precision floating point approximations to decrease the running time 
of the algorithm. Since such an approximation will introduce error, code was included 
to detect and attempt to correct errors resulting from rounding in double-precision 
calculations. In addition, while not changing the operation of the LLL algorithm, 
we performed a change of variables to reduce the rounding error which is necessarily 
introduced. LLL as defined in [27] uses two arrays of variables to hold information 
about the lattice being reduced: 

_ (bi.bj) 

(b,, b,) ’ 

^■ = l|bril = (bf,bf). 
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(Note that the bf are orthogonalized versions of the bi in LLL, not the basis vectors 
of the dual lattice as in Seysen’s reduction algorithm.) Our version of LLL used three 
arrays of variables to maintain the same information: 

(bi.bi) llbjll „ 

miJ ~ INI • ||b;|| M"-’ 

Ci = iibm = jBi 

bb I = (bi,bi) 


It is not difficult to transform the LLL relations between values of fiij and Bi into 
equations involving and bb t . The advantage gained is that: 

(bi, bj) 

llbill-||bl|| ’ 

||bi|| • HbjH • cos a 

INI-||bj|| ’ 

= cos a , 

where a is the angle between the vectors bj and bj. In the original version of LLL, 

_ INI 

A 4 *,.; — ||bt|| 

and thus could take on a wide range of values, depending on the ratio of the lengths 
of bi and bj. Limiting m,-j 6 [—1,1] and keeping the length information separate 
in Ci and bbi reduces the error introduced by double-precision arithmetic operations. 
(Notice that bbi values are integers and can be regenerated at will from the integer 
basis vectors.) Of course, with sufficient multiprecision representation for all variables 
there is no need for this transformation. It is simply an implementation-specific 
modification which allowed us to use a fast, double-precision version of LLL to reduce 
the lattices which arose while solving 106-element subset sum problems. 

Tables 3.1 through 3.3 show the results of experiments performed on random 
subset sum problems with n < 106. For n < 74 twenty subset sum problems were 



3.5. EMPIRICAL TESTS USING ALGORITHM SL 


77 


Table 3.1: Test Results Using Algorithm SL for 42 < n < 58 


n 

b 

d = njb 

RM S 

RMiaAj 

TL 

SstSr) 

% Solved 

Average Time 

20 trials for each length b for each n 

42 

42 

1.000 

36907 

420357 

21 

20 (19) 

1.00 

98.39 

42 

44 

0.955 

38955 

276031 

20 

20 (20) 

1.00 

94.67 

42 

46 

0.913 

40264 

362848 

21 

20 (19) 

1.00 

99.88 

42 

48 

0.875 

42288 

234513 

20 

20 (20) 

1.00 

97.65 

42 

50 

0.840 

43826 

222350 

20 

20 (20) 

1.00 

99.67 

42 

52 

0.808 

45144 

178153 

20 

20 (20) 

1.00 

100.29 

42 

54 

0.778 

46405 

209717 

20 

20 (20) 

1.00 

102.56 

50 

50 

1.000 

55122 

3206204 

42 

20 ( 9) 

1.00 

283.93 

50 

52 

0.962 

56243 

2672637 

36 

19 (10) 

0.95 

262.24 

50 

54 

0.926 

58249 

2112480 

31 

20 (14) 

1.00 

243.22 

50 

56 

0.893 

59721 

1596051 

26 

20 (16) 

1.00 

224.80 

50 

58 

0.862 

61594 

1098460 

21 

20 (19) 

1.00 

208.40 

50 

60 

0.833 

63604 

1220857 

22 

20 (18) 

1.00 

217.68 

50 

62 

0.806 

65563 

978497 

21 

20 (19) 

1.00 

210.71 

58 

58 

1.000 

73587 

15664277 

83 

7(1) 

0.35 

808.68 

58 

60 

0.967 

77371 

13477268 

70 

10(5) 

0.50 

744.32 

58 

63 

0.921 

80574 

13869151 

74 

13(2) 

0.65 

837.96 

58 

66 

0.879 

82284 

12949088 

67 

12(5) 

0.60 

803.07 

58 

69 

0.841 

87569 

7581281 

43 

19 (10) 

0.95 

594.55 

58 

72 

0.806 

90773 

6234735 

38 

20 (10) 

1.00 

548.75 

58 

75 

0.773 

93896 

4408841 

29 

20 (14) 

1.00 

480.70 

58 

78 

0.744 

98289 

3411933 

23 

20 (18) 

1.00 

450.66 

58 

81 

0.716 

101104 

3572104 

25 

20 (15) 

1.00 

463.48 

58 

84 

0.690 

105325 

2090306 

20 

20 (20) 

1.00 

412.65 

58 


0.667 

107711 

1792670 

20 

20 (20) 

1.00 

406.44 


1.00 


406.44 
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attempted per b value. Only ten such problems were attempted per b value for 
82 < n < 106. For each problem, \n elements were chosen from among the weights 
d,... , a n to be in the subset. During the LLL phase of the algorithm, parameter 
d = 5 and x 2 = 8. Six distinct values of y were used for n < 90; nine were used for 
n > 90. 

All tests were run on Silicon Graphics 4D-220 workstations with four or eight 
MIPS Computers, Inc., R3000 processors per workstation. Each workstation was 
equipped with (32Mbytes- # processors) of main memory. As Algorithm SL requires 
significantly less than 32Mbytes of memory per process, all processors in a given 
workstation could be used simultaneously to work on different subset sum problems. 
The majority of the R3000 processors ran at 33MHz; the remainder operated with a 
clock frequency of 25MHz. All running times reported in Tables 3.1 through 3.3 have 
been adjusted to reflect the running time on a 33MHz processor. 

The columns in Tables 3.1 through 3.3 show the value of the following variables 
for each (n, b) pair: 

• n: The number of elements in the set of weights (d>... ,a n ). Exactly one-half 
of these elements were chosen to form the subset sum s. 

• b: The number of bits in the binary representation of the ads. Each a; was 
generated randomly by choosing b random 0 — 1 variables and concatenating 
the bits. 

• d: The density of this class of subset sum problems, d — n/b. 

• RMs: The number of row moves performed by the Seysen phase of Algorithm 
SL on all trials. 

• RMlll'- The number of row moves perform by the LLL phase of Algorithm SL 
on all trials. 
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• TL: The total number of lattices reduced by LLL during all attempted trials, 
where we consider each randomization of the input lattice to be a new lattice. 
Thus, for 20 trials, 20 < TL < 20tti = 100. 

• S 5 : The number of subset sum problems successfully solved with tti = 5. 

• Si: The number of subset sum problems successfully solved with 7 Ti = 1. 
This number is useful when comparing the results of Algorithm SL to the 
Radziszowski-Kreher method. 

• % Solved: The success rate for 7 Ti = 5 measured as a fraction of the number of 
attempted problems. 

• Average Time: The average amount of CPU time (in seconds) required to run 
Algorithm SL on a single trial, adjusted to reflect a 33MHz R3000 processor. 

It should be noted that the running times for Algorithm SL, while an improvement 
over the Radziszowski-Kreher method, could be improved. The code for SL was not 
tuned to the Silicon Graphics machine (although Bailey’s multiprecision code was 
written to work with SGI workstations). Careful recoding and tuning of crucial sub¬ 
routines could probably yield significant improvements without changing the overall 
operation of the algorithm. 

A quick look at the results summarized in Table 3.1 shows that Algorithm SL 
greatly improves upon the performance of previous methods for small values of n 
(n < 58). In [33] Radziszowski and Kreher were able to solve all attempted subset 
sum problems for (n,d) = (42,68), but only managed to solve 30% of the (42,54) 
problems. Compare this result to SL, which was able to solve all attempted subset 
sum problems with n = 42 and density d < 1.0. In fact, SL appears to work as well 
as a lattice oracle for all n < 50. The algorithm first shows signs of degradation at 
n = 58, where a 100% success rate was not reached until b =72. A true lattice oracle 
(using the Euclidean norm) should have been able to achieve 100% success at b = 63, 
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where the density d — 0.926 is less than the critical 0.9408 bound. However, SL was 
still able to solve over half of the attempted subset sum problems at d = 0.962; in 
[33] the 0.50 success rate was not reached until d « 0.56. 

Columns four and five of the Table 3.1 show the number of row moves performed by 
the Seysen and LLL phases of SL. These numbers are a good first-order approximation 
of the amount of work performed by each phase, although it should be pointed out that 
there is no direct correspondence between the absolute magnitudes of the numbers. 
For fixed n, as b increases the Seysen phase performs more reduction steps on the 
lattice basis, and less work is required by the LLL phase to find the solution to 
the subset sum problem. As the density of the system decreases, Seysen’s reduction 
algorithm comes closer to finding the desired solution vector. Indeed, in ten test 
cases with ( n,d ) = (24,24), the LLL phase was able to find solutions using only 
the y = 0.2578125 and y = 0.625 reductions. For higher densities of subset sum 
problems, the Seysen phase of SL reached a local minimum earlier, and LLL had a 
greater reduction to perform to find the solution. Note that once a Seysen reduction 
reached a local minimum it was considered finished; no attempt was made to “kick” 
or “heat” the reduced basis and then reapply Seysen. Application of one or more of 
the techniques described in Section 2.3 might yield better results, and is certainly a 
topic which should be investigated in the future. 

Although Algorithm SL works exceptionally well for small values of n, its per¬ 
formance degrades quickly with respect to n once n > 60. Table 3.2 shows the 
performance of SL for n in the range [66,82]. At n — 66, SL did not obtain a 100% 
success rate until the density dropped to d = 0.66, although the fact that SL solved 
19 out of 20 attempts at b = 92 suggests that incrementing 7Ti and/or 7r 2 might yield 
a 100% success rate at d = 0.717. If we consider only successes on the first attempt, 
SL reached the 100% level at b = 112, compared to b = 144 for Radziszowski and 
Kreher. (Note too that the method in [33] uses -K\ = 1,9 < 7t 2 < 15, whereas in these 
experiments 7r 2 = 8 for SL.) Similar improvements may be seen for the cases where 
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Table 3.2: Test Results Using Algorithm SL for 66 < n < 82 


n 

b 

d = n/b 

RM S 

.RMlll 

TL 

S 5 (Si) 

% Solved 

Average Time 

20 trials for each length b for each n 

66 

76 

0.868 

109671 

36111021 

90 

5(0) 

0.25 

1774.39 

66 

80 

0.825 

116604 

33233636 

84 

7(1) 

0.35 

1670.55 

66 

84 

0.786 

120172 

23212075 

60 

16(5) 

0.80 

1465.55 

66 

88 

0.750 

128015 

24791753 

63 

17(2) 

0.85 

1560.71 

66 

92 

0.717 

134563 

17049844 

45 

19(9) 

0.95 

1249.72 

66 

96 

0.688 

139469 

14280958 

40 

19 (10) 

0.95 

1140.53 

66 

100 

0.660 

142870 

10472664 

30 

20 (11) 

1.00 

993.57 

66 

104 

i 

0.635 

148379 

6026607 

21 

20 (19) 

1.00 

819.64 

66 

108 

0.611 

155150 

5688732 

22 

20 (18) 

1.00 

825.37 

66 

112 

0.589 

162192 

4621556 

20 

20 (20) 

1.00 

803.30 

74 

106 

0.698 

175993 

67316342 

85 

8(2) 

0.40 

3416.51 

74 

112 

0.661 

183540 

55347194 

71 

12(4) 

0.60 

3095.85 

74 

118 

0.627 

189322 

51147665 

66 

13(5) 

0.65 

2850.97 

74 

124 

0.597 

202651 

23245962 

33 

19 (13) 

0.95 

1852.46 

74 

130 

0.569 

211750 

17077901 

25 

20 (16) 

1.00 

1600.57 

10 trials for each length b for each n 

82 

134 

0.612 

122353 

74911469 

50 

0(0) 

0.00 

8499.82 

82 

140 

0.586 

126724 

54569482 

37 

5(1) 

0.50 

6144.08 

82 

146 

0.562 

132948 

63474238 

42 

5(0) 

0.50 

6315.46 

82 

152 

0.539 

139374 

22873751 

17 

10(5) 

1.00 

3375.16 

82 

158 

0.519 

144353 

18199583 

14 

10(7) 

1.00 

3154.05 
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n — 74 and n = 82. Although at n = 82 Algorithm SL did not attain a success rate 
of 1.00 until the density had decreased to d = 0.539, similar results were not reported 
in [33] until d had dropped to below 0.39. 

Table 3.3 shows the performance of SL on subset sum problems with n 100. The 
density at which 100% a success rate is reached is steadily and significantly decreasing 
for small increases in n. By the time n = 106, SL is able to solve all attempted 
problems with density « 0.35, but densities above 0.4 appear unobtainable. 


s 1 



Density 


Figure 3-8: Algorithm SL: Si vs. Density 

Figure 3-8 shows the success rate for solving subset sum problems on the first ran¬ 
domization (Si) versus density for all the test cases described in Tables 3.1 through 3.3. 
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Table 3.3: Test Results Using Algorithm SL for 90 < n < 106 


n 

b 

d = n/b 

RM S 

jRMlll 

TL 

S 5 (Si) 

% Solved 

Average Time 

10 trials for each length b for each n 

— 


0.506 

180739 

111553580 

40 

3(2) 

0.30 

14300.87 

90 

186 

0.484 

188560 

88163144 

32 

8(2) 

0.80 

11241.77 

90 

194 

0.464 

196090 

67630320 

25 

8(4) 

0.80 

9362.65 

90 

202 

0.446 

198909 

48109448 

19 

10(5) 

1.00 

7277.74 

98 

220 

0.445 

238948 

177962307 

38 

3(3) 

0.30 

21725.83 

98 

228 

0.430 

249432 

187347238 

39 

6(0) 

0.60 

21841.45 

98 

236 

0.415 

255213 

125193960 

28 

9(1) 

0.90 

17064.47 

98 

244 

0.402 

264404 

74096833 

17 

10(6) 

1.00 

12386.50 

98 

252 

0.389 

277285 

39984465 

11 

10(9) 

1.00 

9142.08 

98 

260 

0.377 

280876 

42340342 

13 

10(7) 

1.00 

8560.10 

106 

270 

0.393 

322902 

265583840 

37 

5(1) 

0.50 

34147.08 

106 

280 

0.379 

335435 

177869835 

26 

9(2) 

0.90 

25219.92 

106 

290 

0.366 

347679 

175847501 

26 

8(4) 

0.80 

26034.16 

■ 


0.353 

357841 

43319759 

10 

10 (10) 

1.00 

12137.60 

106 

310 

0.342 

369215 

49294587 

12 

10(8) 

1.00 

13184.20 

106 

320 

0.331 

375548 

39512954 

11 

10(9) 

1.00 

11652.82 
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Compare the curves depicted in this figure to those for the Lagarias-Odlyzko (Fig¬ 
ure 3-1) and the Radziszowski-Kreher (Figure 3-3) methods. Algorithm SL shows 
significant improvement over these previous methods and extends the “frontier” of 
solvable subset sum problems. Figure 3-9 shows the success rate versus density when 
SL was allowed to perform five attempted reductions (S5). The improvement made 
by Algorithm SL over previous methods may clearly be seen by comparing this figure 
to those shown previously. The shift in the frontier caused by increasing from one 
to five is also apparent. 


S 5 



Density 


Figure 3-9: Algorithm SL: S 5 vs. Density 


The performance of Algorithm SL is a vast improvement over the techniques used 
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in [26, 33]. Not only is the success rate higher for Algorithm SL, but its improved 
running time allows larger-siged subset sum pr oble m s to be effectively a tta c hed . Nev¬ 
ertheless, for n > 106 Algorithm SL is still only able to attack subset sum problems 
with densities below 0.4. Subset sum problems with larger n and higher densities are 
unaolvable from a practical point of view. 








Chapter 4 


Conclusions 


We have shown in this thesis that Seysen’s lattice basis reduction algorithm performs 
much better than other currently available techniques in a limited number of cir¬ 
cumstances. In particular, we have demonstrated that Seysen’s algorithm may be 
combined with the LLL algorithm to solve subset sum problems. As a general lattice 
basis reduction tool, however, Seysen’s algorithm leaves much to be desired. The 
lack of theoretical bounds on the running time and performance of the algorithm 
is discouraging. Empirical tests performed using Seysen’s algorithm also highlight 
weaknesses with using this method to solve general lattice reduction problems. As 
the first part of this thesis demonstrated, the performance of Seysen’s algorithm on 
randomly generated “dense” lattice bases degrades quickly as the dimension of the 
lattice increases above a certain critical bound. There are numerous local minima of 
the S'(A) function for these lattices, and once Seysen’s algorithm encounters one it is 
unable to escape (without some external influence acting upon the lattice basis). 

Seysen’s algorithm should not be immediately discounted, however. In certain 
specific cases the algorithm performed quite well. In general Seysen’s algorithm per¬ 
forms only a fraction of the row moves required by the LLL algorithm to reduce a 
lattice basis. For certain lattices, such as the Be lattices of Section 2.2.3, the LL- 
L algorithm was unable to perform any row reductions, whereas Seysen’s algorithm 
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enjoyed great success. Even for dense lattices, if the dimension of the lattice is rel¬ 
atively small, Seysen’s algorithm obtains the same results as LLL in only a fraction 
of the time. Applications which utilize basis reduction in few dimensions are good 
candidates for Seysen’s method. 

The experiments performed on lattices derived from subset sum problems highlight 
one of the main advantages of Seysen’s technique: its ability to very quickly perform 
significant reductions on a lattice basis. Using Seysen’s algorithm as a first reduction 
stage allowed us to convert lattice bases involving large, multiprecision values into 
other bases with vector coefficients which could be represented in single- or double¬ 
precision. The LLL algorithm could then be run without having to use multiprecision 
arithmetic, which greatly improved its execution time. For larger lattices (with n « 
100 ) one must be aware of possible problems due to rounding and truncation errors, 
but these difficulties can be overcome. 

Although the primary goal of this thesis was to investigate applications of Seysen’s 
basis reduction algorithm, one should not overlook the other techniques which were 
used (in addition to Seysen’s algorithm) to solve subset sum problems. In particular, 
the theoretical improvements to the Lagarias-Odlyzko attack in [12] may be directly 
incorporated into practical methods of solving subset sum problems. Also, the use of 
multiple values of the y parameter in the LLL algorithm significantly reduced the total 
number of row moves LLL performed. This modification does not appear to decrease 
the reduction performance of LLL in any way over LLL with constant y ~ 0.99, 
although it does significantly reduce the running time of the algorithm. We would 
suggest using multiple, increasing values of y in the future whenever LLL-reduction 
is performed. 

We have demonstrated that Seysen’s algorithm is a good basis reduction technique 
for certain types of lattices, outperforming the basic LLL algorithm in terms of the 
number of row operations required. Furthermore, we have seen how Seysen’s algo¬ 
rithm may be combined with variants of the LLL algorithm and heuristic methods 
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to successfully attack many subset sum problems with n < 100. Yet these specific 
instances are but a fraction of the type of lattice reduction problems which arise. It 
is therefore natural to consider what other types of lattice reduction problems are 
suitable for attack by Seysen’s algorithm. We conclude this chapter, and the the¬ 
sis itself, with some remarks on other Seysen-suitable lattices, and also suggestions 
for modifying Algorithm SL to solve subset sum problems of larger dimensions and 
higher densities. 


4.1 Candidate Lattices for Seysen Reduction 


The number of “classes” of lattice bases which may be successfully reduced using only 
Seysen’s basis reduction algorithm appears to be quite small. It was shown in Chap¬ 
ter 2 that Seysen’s algorithm, which was designed for finding simultaneously good 
reductions of a lattice basis and its dual, indeed works well in such cases. Further¬ 
more, if the goal of a lattice reduction is to minimize some cost function or measure 
of the lattice (and perhaps its dual), an appropriately modified version of Seysen’s 
algorithm incorporating the cost function will likely perform much better than an 
LLL-type algorithm. While there are empirical reasons which suggest the use of the 
5(A) measure as a particular cost function, Seysen himself mentioned in [38] that 
other function had been used in the algorithm in place of 5(A) with about the same 
degree of success. 

Replacing 5(A) with another cost function may involve some difficulty. In par¬ 
ticular, it may not be possible to solve for the coefficients A,- in closed form using an 
alternate cost function. If some range of acceptable values can be determined, then 
a search considering all integers in the range of interest may be feasible. Such an 
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approach was used in experiments with the cost function 

i=l Y 
n i 

= £ 


U IN 

For a proposed row move bj «— bj + Ab, it is possible to compute bounds on the 
range of possible A values. For sufficiently small ranges A 6 [Ax, A 2 ] we compute the 
change in 5' / (>l) after applying the row move with A taking on every integer value 
in the range. Then (assuming a greedy approach) we choose the value for A which 
maximized the decrease in .^'(A). Thus, although it may not be possible to solve for 
the best choice for A for some proposed cost function, numerical methods may allow 
consideration of the metric anyway. 

The class of lattice bases which may be successfully reduced by using only Seysen’s 
algorithm appears quite limited, assuming that we ignore differences in execution 
time. The randomly generated lattice bases with unit determinant from Chapter 2 
exemplify this point: Seysen’s algorithm began to perform poorly for n > 35, whereas 
the LLL algorithm continued to correctly reduce lattices to the n-dimensional cubic 
lattice for all tested cases. As a stand-alone technique, Seysen’s algorithm may not be 
considered that useful; it works quite well for some lattices of low dimension, but it 
tends to stop early for larger-dimensioned lattice bases. To reduce these lattice bases 
some combination of Seysen’s algorithm and other reduction methods is probably 
required. 


4.2 Modifying Algorithm SL 

The combination of the Seysen and LLL reduction algorithms used in Chapter 3 to 
solve subset sum problems is a significant improvement over previously tried tech¬ 
niques. There is still plenty of room for improvement. The performance of Algorithm 
SL declines steadily over the range 50 < n < 100. For subset sum problems with 
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n m 100, Algorithm SL had difficulty solving instances with density d « 0.4, well 
below the 0.9408 ... theoretical density bound. 

There are a number of ways in which Algorithm SL could be extended in order to 
attack subset sum problems with larger n, d values. For example, any of the techniques 
mentioned in Section 2.3 above could be applied to the Seysen phases of the algorithm. 
We can clearly imagine that another measure function S'(A) might perform better for 
sparse, subset-sum generated lattices when compared to the performance of .S'(A) = 

EJU INHIbril 2 . 

Another method of extending Algorithm SL would be to increase the values of 
the 7Ti and 7T2 parameters. We have seen that increasing Xi (the number of different 
orderings of the initial LLL lattice basis vectors) from one to five yields markedly 
better performance. Along similar lines, one could also allow more LLL- Weight- 
Reduction-Sort loops to occur for each randomization. Recall that for Algorithm SL, 
7T2, the maximum number of loops, was set to eight. After eight iterations without 
finding the desired solution vector in the lattice basis, Algorithm SL would “give up” 
and select a new randomization of the LLL-Phase input lattice basis. The choice 
to set 7T2 = 8 was made arbitrarily; Radziszowski and Kreher used 7r2 = 15 in their 
experiments. The majority of the work performed during the LLL-Phase of Algorithm 
SL occurs during the first LLL-Loop on each of the tti iterations; the second through 
7T2 h LLL-Loops perform only a fraction of the number of row moves performed by the 
first loop. In our experiments, subsequent LLL-Loops tended to perform only ^ the 
row operations performed in the initial LLL- Weight-Reduction-Sort iteration. This 
means that i r 2 could probably be increased to around 20 and the overall running time 
of the LLL-Phase would probably double. 

Increasing the number of LLL-Loops appears to be an inexpensive way to allow 
Algorithm SL to perform more reduction operations on a lattice basis. One must 
consider, however, how much an increase in 7r 2 will improve the rate at which Al¬ 
gorithm SL solves subset sum problems. We know from [33] that Radziszowski and 
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Kreher were able to improve upon the Lagarias-Odlyzko algorithm by significantly 
increasing 7r 2 . Furthermore, in some very recent experiments Euchner and Schnorr 
(personal communication) were able to obtain performance close to that of Algorithm 
SL with tti = 1 and 7 t 2 = 16. The Euchner-Schnorr reduction algorithm utilizes some 
other heuristics not included in Algorithm SL but uses only the LLL algorithm as 
the main reduction technique. To date, Euchner and Schnorr have only reported on 
experiments involving subset sum problems of size n < 66; the performance of their 
technique on lattice bases in higher dimensions is unknown. Even their preliminary 
results, though, suggest that it might be worthwhile to increase the number of reduc¬ 
tion stages in Algorithm SL, even if it is necessary to reduce in order to maintain 
comparable running times. Also, as their methods uses only the LLL algorithm, a 
combination of their techniques and Algorithm SL is certainly possible. 

Another factor to consider in Algorithm SL is the structure of our version of LLL. 
Recall that instead of running LLL with y < 1.0 we used a fixed number of LLL 
stages with varying y values. The first stage used a value of y slightly greater than 1, 
and subsequent stages used larger and larger values of y. For the first stage of the i r 2 
reduction stages, the multiple-value version of LLL significantly reduced the overall 
running time of the algorithm. However, later reduction stages were generally unable 
to benefit from this structure; in fact, in many cases the lattice was reduced further 
only by the stages with y > 0.875. Future implementations might consider removing 
the LLL stages with small y values for the second through (7r 2 ) th reduction loops. 

In Section 2.3.3 we discussed how Hadamard matrices could be used to randomly 
permute a lattice basis which was locally Seysen-reduced in the hope that reapplying 
Seysen’s algorithm would yield a better reduced lattice basis. This same technique 
could be applied to both the Seysen and LLL phases of Algorithm SL and might yield 
significant improvements. Hadamard matrices could be used to permute a Seysen- 
reduced lattice basis whose Seysen measure is a local minimum of the S(A ) function. 
Such a permutation might then permit Seysen’s algorithm to reach another lattice 
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basis with smaller S'(A) value, which in turn could increase the ratio of work performed 
by the Seysen phase to that performed by the LLL phase. Similarly, recall that 
under the current scheme once a lattice basis has been LLL-reduced 7T2 times without 
yielding a solution vector, that basis is forgotten and a new iteration is started using 
a random permutation of the output of the Seysen phase. Instead of “throwing away” 
the LLL-reduced basis, which has shorter basis vectors than the Seysen-reduced basis, 
a Hadamard permutation could be applied to the output of the tt 2 stages. This 
permuted basis would then be used as the input basis to the next big iteration. 
Assuming that the Hadamard permutation method sufficiently “scrambles” the lattice 
basis, we could thus avoid the overhead of having to reduce the output of the Seysen 
phase more than once. 

Finally, it is possible to modify Algorithm SL to take into account additional 
information related to the creation of the specific subset sum problem. In [12] it is 
shown how to tailor the input lattice basis for a subset sum problem if it is known 
that the number of elements in the desired subset of weights is bounded by a fraction 
of n. That is, if it is known that 

n 

Yj e i = P n i 

t = l 

then the lattice basis can be modified so that a solution vector exists of length 
\Jnf3{ 1 — /?). (In the worst case, j3 = \ and the solution vector is the familiar e' vector 
with He'll = For instances of the general subset sum problem no information is 

known concerning Some knapsack cryptosystem, such as the Chor-Rivest sys¬ 

tem [11], do use subsets with relatively few weights. When attacking such systems, 
Algorithm SL should be modified to use the tailored lattice basis described in [12]. 
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